package bs.bs2d.fea;
import bs.bs2d.geom.JTSUtils;
import bs.bs2d.gui.UserMessage;
import bs.bs2d.io.CalculixIO;
import bs.bs2d.io.MeshReader;
import com.vividsolutions.jts.geom.Coordinate;
import com.vividsolutions.jts.geom.LineString;
import com.vividsolutions.jts.geom.MultiLineString;
import com.vividsolutions.jts.geom.MultiPolygon;
import com.vividsolutions.jts.geom.Point;
import java.awt.Shape;
import java.awt.geom.Rectangle2D;
import java.io.File;
import java.util.List;
import javax.swing.filechooser.FileFilter;
import javax.swing.filechooser.FileNameExtensionFilter;
import org.apache.commons.io.FileUtils;
import org.apache.commons.io.FilenameUtils;
/**
* Class to contain and manipulate result data of a finite element analysis
* (displacement and stresses) on a 2D tri-mesh.
*
* @author Djen
*/
public class FEAResults {
private static FileFilter fileFilter = null;
public enum Scalar{
VON_MISES_STRESS,
X_DISPLACEMENT,
Y_DISPLACEMENT,
TOTAL_DISPLACEMENT,
LOCAL_STRAIN,
MAX_PRINCIPAL_STRAIN,
MAX_ABSOLUTE_STRAIN,
DISTANCE_FROM_PERIMETER
}
private TriMesh2D mesh;
private float[] xDisplacement;
private float[] yDisplacement;
private float[][] elementStress;
public static FEAResults getResultsFromFile(TriMesh2D mesh, File f){
if(f.exists()){
float[][] disp = CalculixIO.readDisplacements(f);
float[][] stress = CalculixIO.readStress(f);
return new FEAResults(mesh, disp, stress);
}
UserMessage.showError("Result File Missing", "The specified file does not exist!");
return null;
}
/**
* Retrieves FEA results from a *.dat file. This method assumes that the
* associated mesh can be found in a *.msh file with the same name and in
* the same location as the *.dat file.
* @param dat the *.dat file to read results from
* @return the FEA results read from dat
*/
public static FEAResults getResultsFromFile(File dat){
File msh = FileUtils.getFile(FilenameUtils.removeExtension(dat.getPath()) + ".msh");
if (!msh.exists()) {
UserMessage.showError("Mesh file missing!", "The results file is missing the associated mesh.");
return null;
}
try {
TriMesh2D mesh = MeshReader.readMesh(msh);
return getResultsFromFile(mesh, dat);
} catch (Exception ex) {
UserMessage.showError("Unable to load results!", "The result file could not be loaded.");
ex.printStackTrace();
return null;
}
}
/**
* Creates a new FEAResult instance.
* @param mesh the mesh to which the analysis values respond
* @param displacement the displacement in pairs of x and y values
* @param elementStress the stress values for each element (sxx, syy, szz, sxy,sxz,syz)
*/
public FEAResults(TriMesh2D mesh, float[][] displacement, float[][] elementStress) {
if(displacement == null || elementStress == null ||
displacement.length != mesh.getNodeCount() ||
elementStress.length != mesh.getElementCount()){
throw new IllegalArgumentException("Result data does not fit mesh. Check result file!");
}
this.mesh = mesh;
this.elementStress = elementStress;
splitDisplacements(displacement);
}
private void splitDisplacements(float[][] disp){
xDisplacement = new float[disp.length];
yDisplacement = new float[disp.length];
for (int i = 0; i < disp.length; i++) {
xDisplacement[i] = disp[i][0];
yDisplacement[i] = disp[i][1];
}
}
/**
* @return the original mesh the FEA was performed on
*/
public TriMesh2D getMesh() {
return mesh;
}
/**
* Returns a deformed copy of the mesh. Each nodes is moved according to
* the displacement result values, scaled by the given scale factor.
* @param scale a factor by wich to scale the displacement
* @return a deformed mesh
*/
public TriMesh2D getDeformedMesh(float scale){
TriMesh2D dm = mesh.copy();
List<Node2D> nodes = dm.getNodes();
for(Node2D n : nodes){
int i = n.getIndex();
n.move(scale * xDisplacement[i], scale * yDisplacement[i]);
}
return dm;
}
/**
* Computes a displacement scale factor based on the size of the original
* mesh and a given relative deformation factor. Can be used to
* automatically scale the deformation to visible displacements.
* @param relDef a relative deformation factor. Would typically be something
* around 0.1
* @return a scale factor to use with the getDeformedMesh method
*/
public float getMeshDeformationScale(float relDef){
float maxX = 0;
float maxY = 0;
for (int i = 0; i < xDisplacement.length; i++) {
float x = Math.abs(xDisplacement[i]);
float y = Math.abs(yDisplacement[i]);
if(x > maxX)
maxX = x;
if(y > maxY)
maxY = y;
}
Rectangle2D bounds = mesh.getBounds();
float sx = (float)bounds.getWidth() * relDef / maxX;
float sy = (float)bounds.getHeight() * relDef / maxY;
return Math.min(sx, sy);
}
public float[] getScalar(Scalar scalar){
switch(scalar){
case VON_MISES_STRESS:
return interpolateToNodes(mesh, getVonMisesStress(elementStress));
case X_DISPLACEMENT:
return xDisplacement;
case Y_DISPLACEMENT:
return yDisplacement;
case TOTAL_DISPLACEMENT:
return getTotalDisplacement(xDisplacement, yDisplacement);
case LOCAL_STRAIN:
// return computeNodeStrain(mesh, getDeformedMesh(1));
return interpolateToNodes(mesh, computeElementStrain(mesh, getDeformedMesh(1)));
case MAX_PRINCIPAL_STRAIN:
return interpolateToNodes(mesh, getMaxPrincipalStrain(elementStress));
case MAX_ABSOLUTE_STRAIN:
return interpolateToNodes(mesh, getMaxAbsoluteStrain(elementStress));
case DISTANCE_FROM_PERIMETER:
return getDistanceFromPerimeter(mesh);
default:
throw new IllegalArgumentException("Unknown scalar!");
}
}
private static float[] getTotalDisplacement(float[] xDisp, float[] yDisp){
float[] td = new float[xDisp.length];
for (int i = 0; i < td.length; i++) {
td[i] = (float)Math.sqrt(xDisp[i] * xDisp[i] + yDisp[i] * yDisp[i]);
}
return td;
}
private static float[] getDistanceFromPerimeter(TriMesh2D mesh){
float[] d = new float[mesh.getNodeCount()];
Shape shape = mesh.createOutlinePath();
MultiPolygon mp = JTSUtils.createMultiPolygon(shape);
LineString[] ls = JTSUtils.getLineStrings(mp);
MultiLineString mls = JTSUtils.getFactory().createMultiLineString(ls);
for (int i = 0; i < d.length; i++) {
Node2D n = mesh.getNode(i);
Coordinate c = new Coordinate(n.getX(), n.getY());
Point p = JTSUtils.getFactory().createPoint(c);
d[i] = (float) p.distance(mls);
}
return d;
}
private static float[] computeElementStrain(TriMesh2D mesh, TriMesh2D defMesh){
float[] strain = new float[mesh.getElementCount()];
for (int i = 0; i < mesh.getElementCount(); i++) {
float area = mesh.getElement(i).getArea();
float defArea = defMesh.getElement(i).getArea();
strain[i] = defArea / area - 1;
}
return strain;
}
@Deprecated
private static float[] computeEdgeStrain(TriMesh2D mesh, TriMesh2D defMesh){
float[] strain = new float[mesh.getEdgeCount()];
for (int i = 0; i < mesh.getEdgeCount(); i++) {
Edge edge = mesh.getEdge(i);
Edge defEdge = defMesh.getEdge(i);
strain[i] = defEdge.getLength() / edge.getLength() - 1;
}
return strain;
}
@Deprecated
private static float[] computeNodeStrain(TriMesh2D mesh, TriMesh2D defMesh){
float[] eStrain = computeEdgeStrain(mesh, defMesh);
float[] nStrain = new float[mesh.getNodeCount()];
for (int i = 0; i < mesh.getNodeCount(); i++) {
Node2D n = mesh.getNode(i);
float val = 0;
for (int j = 0; j < n.getEdgeCount(); j++) {
val += eStrain[n.getEdge(j).getIndex()];
}
nStrain[i] = val / n.getEdgeCount();
}
return nStrain;
}
/**
* Computes von mises stress from principal and shear stress for each element.
* @param stress the stress componenets for each element
* @return the von mises stress values for each element
*/
private static float[] getVonMisesStress(float[][] stress){
float[] mises = new float[stress.length];
for (int i = 0; i < mises.length; i++) {
float[] d = stress[i];
// calculate von mises stress
mises[i] = (float) Math.sqrt(d[0] * d[0] + d[1] * d[1]
- d[0] * d[1] + d[2] * d[2] * 3);
}
return mises;
}
private static float[] getMaxPrincipalStrain(float[][] stress){
float u = 0.3f; //poisson ratio
float[] mps = new float[stress.length];
for(int i = 0; i < mps.length; i++){
float[] d = stress[i];
// calculate maximum principal strain
mps[i] = Math.max(d[0] - u*d[1], d[1] - u*d[0]);
}
return mps;
}
private static float[] getMaxAbsoluteStrain(float[][] stress){
float u = 0.3f; //poisson ratio
float[] mas = new float[stress.length];
for(int i = 0; i < mas.length; i++){
float[] d = stress[i];
// calculate maximum principal strain
float dx = d[0] - u*d[1];
float dy = d[1] - u*d[0];
// compare absolutes
if(Math.abs(dx) > Math.abs(dy))
mas[i] = dx;
else
mas[i] = dy;
}
return mas;
}
/**
* Interpolates node data from element data by computing the average of the
* values of all adjacent elements of each node.
*/
private static float[] interpolateToNodes(TriMesh2D mesh, float[] elementScalar){
float[] nodeScalar = new float[mesh.getNodeCount()];
// interpolate stress in nodes
// --> averaging stress in adjacent elements
for(Node2D n : mesh.getNodes()){
float avrg = 0;
for(TriElement e : n.getElements()){
avrg += elementScalar[e.getIndex()];
}
avrg /= n.getElementCount();
nodeScalar[n.getIndex()] = avrg;
}
return nodeScalar;
}
public static FileFilter getFileFilter() {
if(fileFilter == null)
fileFilter = new FileNameExtensionFilter("Calculix output file (*.dat)", "dat");
return fileFilter;
}
}