Skip to content

Commit d5df2c1

Browse files
committed
Moved some tools to separate files to make them available across more modifiers
1 parent 28353ae commit d5df2c1

File tree

2 files changed

+217
-0
lines changed

2 files changed

+217
-0
lines changed

FEBioStudio/MeasureTools.cpp

Lines changed: 183 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,183 @@
1+
/*This file is part of the FEBio Studio source code and is licensed under the MIT license
2+
listed below.
3+
4+
See Copyright-FEBio-Studio.txt for details.
5+
6+
Copyright (c) 2025 University of Utah, The Trustees of Columbia University in
7+
the City of New York, and others.
8+
9+
Permission is hereby granted, free of charge, to any person obtaining a copy
10+
of this software and associated documentation files (the "Software"), to deal
11+
in the Software without restriction, including without limitation the rights
12+
to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
13+
copies of the Software, and to permit persons to whom the Software is
14+
furnished to do so, subject to the following conditions:
15+
16+
The above copyright notice and this permission notice shall be included in all
17+
copies or substantial portions of the Software.
18+
19+
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
20+
IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
21+
FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
22+
AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
23+
LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
24+
OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
25+
SOFTWARE.*/
26+
27+
// This calculates the moment of inertia, but only approximately!
28+
// Requires a fine mesh for accurate results.
29+
30+
#include "stdafx.h"
31+
#include <MeshLib/FSMeshBuilder.h>
32+
#include <MeshLib/MeshTools.h>
33+
#include <FECore/mat3d.h>
34+
#include <FECore/vec3d.h>
35+
#include "MeasureTools.h"
36+
#include <GeomLib/GObject.h>
37+
38+
// This algorithm calculates the center of mass approximately,
39+
// by weighing the element centers by the element volumes.
40+
// TODO: I should look into bringing the numerical integration tools
41+
// from FEBio to FBS so this can be calculated accurately.
42+
vec3d CalculateCOM(FSMesh& mesh)
43+
{
44+
vec3d rc(0, 0, 0);
45+
double Vtotal = 0.0;
46+
47+
int selectedElements = mesh.CountSelectedElements();
48+
49+
// loop over all elements
50+
int NE = mesh.Elements();
51+
for (int iel = 0; iel < NE; ++iel)
52+
{
53+
FSElement& el = mesh.Element(iel);
54+
if ((selectedElements == 0) || el.IsSelected())
55+
{
56+
// Volume is calculated in global frame
57+
double V0 = mesh.ElementVolume(el);
58+
59+
// center of mass is calculated in local frame
60+
vec3d c_local = mesh.ElementCenter(el);
61+
vec3d c_global = mesh.LocalToGlobal(c_local);
62+
63+
rc += c_global * V0;
64+
Vtotal += V0;
65+
}
66+
}
67+
if (Vtotal != 0.0)
68+
rc /= Vtotal;
69+
70+
const double eps = 1e-12;
71+
double L = rc.Length();
72+
if (L > eps)
73+
{
74+
if (fabs(rc.x) < eps) rc.x = 0;
75+
if (fabs(rc.y) < eps) rc.y = 0;
76+
if (fabs(rc.z) < eps) rc.z = 0;
77+
}
78+
79+
return rc;
80+
}
81+
82+
vec3d CalculateAreaCOM(FSMesh& mesh)
83+
{
84+
int selectedFaces = mesh.CountSelectedFaces();
85+
vec3d c = vec3d(0, 0, 0);
86+
double Atotal = 0.0;
87+
88+
// loop over all elements
89+
for (int i = 0; i < mesh.Faces(); ++i)
90+
{
91+
FSFace& face = mesh.Face(i);
92+
if ((selectedFaces == 0) || face.IsSelected())
93+
{
94+
double Ai = mesh.FaceArea(face);
95+
vec3d ci = mesh.FaceCenter(face);
96+
97+
Atotal += Ai;
98+
c += ci * Ai;
99+
}
100+
}
101+
if (Atotal != 0.0)
102+
c /= Atotal;
103+
104+
return mesh.LocalToGlobal(c);
105+
}
106+
107+
mat3d CalculateMOI(FSMesh& mesh)
108+
{
109+
int selectedElements = mesh.CountSelectedElements();
110+
mat3dd I(1); // identity tensor
111+
112+
mat3d moi; moi.zero();
113+
114+
// calculate the COM
115+
vec3d com = CalculateCOM(mesh);
116+
117+
// loop over all elements
118+
int NE = mesh.Elements();
119+
for (int iel = 0; iel < NE; ++iel)
120+
{
121+
FSElement& el = mesh.Element(iel);
122+
123+
if ((selectedElements == 0) || el.IsSelected()) {
124+
// ElementCenter is calculated in local frame
125+
vec3d c_local = mesh.ElementCenter(el);
126+
vec3d c_global = mesh.LocalToGlobal(c_local);
127+
128+
// Element volume is calculated in global frame
129+
double V0 = mesh.ElementVolume(el);
130+
131+
vec3d r = c_global - com;
132+
mat3d Iij = (r * r) * I - (r & r);
133+
moi += Iij * V0;
134+
}
135+
}
136+
137+
const double eps = 1e-12;
138+
return CleanUp(moi, eps);
139+
}
140+
141+
mat3d CalculateAreaMOI(FSMesh& mesh)
142+
{
143+
int selectedFaces = mesh.CountSelectedFaces();
144+
mat3dd I(1); // identity tensor
145+
146+
mat3d moi; moi.zero();
147+
148+
// calculate the COM
149+
vec3d com = CalculateAreaCOM(mesh);
150+
151+
// loop over all elements
152+
for (int i = 0; i < mesh.Faces(); ++i)
153+
{
154+
FSFace& face = mesh.Face(i);
155+
if ((selectedFaces == 0) || face.IsSelected()) {
156+
double Ai = mesh.FaceArea(face);
157+
// face center is calculated in local frame
158+
vec3d ci = mesh.FaceCenter(face);
159+
160+
vec3d cg = mesh.LocalToGlobal(ci);
161+
162+
vec3d r = cg - com;
163+
mat3d Iij = (r * r) * I - (r & r);
164+
moi += Iij * Ai;
165+
}
166+
}
167+
168+
const double eps = 1e-12;
169+
return CleanUp(moi, eps);
170+
}
171+
172+
mat3d CleanUp(mat3d& m, const double eps)
173+
{
174+
double L = m.norm();
175+
if (L > eps)
176+
{
177+
for (int i = 0; i < 3; ++i)
178+
for (int j = 0; j < 3; ++j)
179+
if (fabs(m[i][j]) < eps) m[i][j] = 0.0;
180+
}
181+
return m;
182+
183+
}

FEBioStudio/MeasureTools.h

Lines changed: 34 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,34 @@
1+
/*This file is part of the FEBio Studio source code and is licensed under the MIT license
2+
listed below.
3+
4+
See Copyright-FEBio-Studio.txt for details.
5+
6+
Copyright (c) 2025 University of Utah, The Trustees of Columbia University in
7+
the City of New York, and others.
8+
9+
Permission is hereby granted, free of charge, to any person obtaining a copy
10+
of this software and associated documentation files (the "Software"), to deal
11+
in the Software without restriction, including without limitation the rights
12+
to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
13+
copies of the Software, and to permit persons to whom the Software is
14+
furnished to do so, subject to the following conditions:
15+
16+
The above copyright notice and this permission notice shall be included in all
17+
copies or substantial portions of the Software.
18+
19+
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
20+
IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
21+
FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
22+
AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
23+
LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
24+
OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
25+
SOFTWARE.*/
26+
27+
#pragma none
28+
29+
vec3d CalculateCOM(FSMesh& mesh);
30+
vec3d CalculateAreaCOM(FSMesh& mesh);
31+
mat3d CalculateMOI(FSMesh& mesh);
32+
mat3d CalculateAreaMOI(FSMesh& mesh);
33+
mat3d CleanUp(mat3d& m, const double eps);
34+

0 commit comments

Comments
 (0)