Skip to content

Commit c797198

Browse files
epsi1onepsi1on
committed
minor update
1 parent f979a34 commit c797198

6 files changed

Lines changed: 181 additions & 7 deletions

File tree

‎BriefFiniteElementNet.TestConsole/Program.cs‎

Lines changed: 3 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -36,13 +36,15 @@ static void Main(string[] args)
3636
{
3737
Console.Title = "BFE tests & temporary codes";
3838

39+
Validation.GithubIssues.Issue96.Run1();
40+
3941
TestTet();
4042
//TestTriangle();
4143
//Validation.GithubIssues.Issue50.Run1();
4244
new Validation.Case_01.Validator().Validate();
4345
//TestGrid();
4446
return;
45-
Validation.GithubIssues.Issue41.Run3();
47+
4648

4749
//TestHingedInternalForce();
4850
//TestBinModel();

‎BriefFiniteElementNet.Validation/BriefFiniteElementNet.Validation.csproj‎

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -63,6 +63,7 @@
6363
<Compile Include="Case_01\Validator.cs" />
6464
<Compile Include="Case_03\Validator.cs" />
6565
<Compile Include="Case_04\Validator.cs" />
66+
<Compile Include="Case_05\Validator.cs" />
6667
<Compile Include="EqualDofValidator.cs" />
6768
<Compile Include="GithubIssues\Issue23.cs" />
6869
<Compile Include="GithubIssues\Issue24.cs" />
@@ -71,6 +72,7 @@
7172
<Compile Include="GithubIssues\Issue41.cs" />
7273
<Compile Include="GithubIssues\Issue48.cs" />
7374
<Compile Include="GithubIssues\Issue50.cs" />
75+
<Compile Include="GithubIssues\Issue96.cs" />
7476
<Compile Include="OpenseesTclGenerator\TetrahedralElementToTcl.cs" />
7577
<Compile Include="OpenseesTclGenerator\TetrahedronToTcl.cs" />
7678
<Compile Include="OpenseesTclGenerator\TriangleElementToTcl.cs" />

‎BriefFiniteElementNet.Validation/Case_03/Validator.cs‎

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -17,7 +17,7 @@ public ValidationResult Validate()
1717

1818
/**/
1919
{
20-
var model = StructureGenerator.Generate3DTetrahedralElementGrid(4, 4, 100);
20+
var model = StructureGenerator.Generate3DTetrahedralElementGrid(3, 3, 70);
2121

2222
var e = 210e9;
2323

@@ -50,15 +50,15 @@ public ValidationResult Validate()
5050
rigid.Nodes.Add(node);
5151
}
5252

53-
//model.MpcElements.Add(rigid);
53+
model.MpcElements.Add(rigid);
5454
model.Trace.Listeners.Add(new ConsoleTraceListener());
5555
model.Solve_MPC();
5656

5757

5858
var delta = f * l * l * l / (3 * e * I);
5959

6060
var t = cnt.FirstOrDefault().GetNodalDisplacement();
61-
61+
6262
var ratio = delta / t.DX;
6363
}
6464
/**/

‎BriefFiniteElementNet.Validation/Case_05/Validator.cs‎

Lines changed: 112 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,112 @@
1+
using System;
2+
using System.Collections.Generic;
3+
using System.Globalization;
4+
using System.Linq;
5+
using System.Text;
6+
using System.Threading.Tasks;
7+
using BriefFiniteElementNet.Common;
8+
using BriefFiniteElementNet.Elements;
9+
using BriefFiniteElementNet.MpcElements;
10+
using HtmlTags;
11+
12+
namespace BriefFiniteElementNet.Validation.Case_05
13+
{
14+
[ValidationCase("Console beam with tetrahedron", typeof(RigidElement_MPC))]
15+
public class Validator : IValidationCase
16+
{
17+
public ValidationResult Validate()
18+
{
19+
var model = StructureGenerator.Generate3DTetrahedralElementGrid(2, 2, 50);
20+
21+
RigidElement_MPC mpc;
22+
23+
{
24+
var e = 210e9;
25+
26+
27+
foreach (var elm in model.Elements)
28+
{
29+
if (elm is TetrahedronElement)
30+
{
31+
var tet = elm as TetrahedronElement;
32+
33+
tet.Material = new Materials.UniformIsotropicMaterial(e, 0.25);
34+
}
35+
}
36+
37+
var dx = model.Nodes.Max(i => i.Location.X) - model.Nodes.Min(i => i.Location.X);
38+
var dy = model.Nodes.Max(i => i.Location.Y) - model.Nodes.Min(i => i.Location.Y);
39+
var dz = model.Nodes.Max(i => i.Location.Z) - model.Nodes.Min(i => i.Location.Z);
40+
41+
var l = dz;// model.Nodes.Max(i => i.Location.Z);
42+
43+
var cnt = model.Nodes.Where(i => i.Location.Z == l);
44+
45+
var f = 1e7;
46+
var I = dy * dx * dx * dx / 12;
47+
var rigid = mpc = new MpcElements.RigidElement_MPC() {UseForAllLoads = true};
48+
49+
foreach (var node in cnt)
50+
{
51+
node.Loads.Add(new NodalLoad(new Force(f / cnt.Count(), 0, 0, 0, 0, 0)));
52+
rigid.Nodes.Add(node);
53+
}
54+
55+
model.MpcElements.Add(rigid);
56+
model.Trace.Listeners.Add(new ConsoleTraceListener());
57+
model.Solve_MPC();
58+
59+
60+
var delta = f * l * l * l / (3 * e * I);
61+
62+
var t = cnt.FirstOrDefault().GetNodalDisplacement();
63+
64+
var ratio = delta / t.DX;
65+
}
66+
67+
68+
{
69+
var span = new HtmlTag("span");
70+
span.Add("p").Text("Validation of rigid element");
71+
span.Add("paragraph").Text("The rigid element connects different nodes through non-deformable elements. ");
72+
span.Add("h3").Text("Model Definition");
73+
74+
span.Add("paragraph").Text(string.Format(CultureInfo.CurrentCulture, "A console beam with tetrahedron, at the end (tip of console) nodes are connected with a rigid element.")).AddClosedTag("br");
75+
76+
span.Add("paragraph").Text("The start of beam is totally fixed in the 3D space.").AddClosedTag("br");
77+
span.Add("paragraph").Text("The end of beam is loaded with a vertical load.").AddClosedTag("br");
78+
span.Add("paragraph").Text("All nodes connected to rigid element should have a rigid body displacement.").AddClosedTag("br");
79+
80+
var mpcNodes = mpc.Nodes.ToArray();
81+
82+
var disps = mpcNodes.Select(i => Tuple.Create(i.Location, i.GetNodalDisplacement()));
83+
84+
85+
span.Add("h3").Text("Validation Result");
86+
87+
if (disps.Any(i => double.IsNaN(i.Item2.DX)|| double.IsNaN(i.Item2.DY)|| double.IsNaN(i.Item2.DZ)|| double.IsNaN(i.Item2.RX)))
88+
{
89+
span.Add("paragraph").Add("h1").Text(string.Format("Validation failed, there are NaNs in the result for nodal displacements"));
90+
}
91+
92+
span.Add("paragraph").Text(string.Format("Validation output for nodal displacements:"));
93+
var n1 = model.Nodes[1].GetNodalDisplacement();
94+
var n2 = model.Nodes[2].GetNodalDisplacement();
95+
var diff = n1 - n2;
96+
//linear beam with loads -> rigid elements introduces rigid body -> Y and Z (translations perpendicular to the beam) are non zero and scaled. Set to zero to avoid confusion with result.
97+
diff.DY = diff.DZ = 0.0;
98+
99+
//span.Add("p").AddClass("bg-info").Text(string.Format("-Max ABSOLUTE Error: {0:e3}", diff));//htmltags cannot encode the delta and teta chars so will use vector length
100+
span.Add("p").AddClass("bg-info")
101+
.Text(string.Format("-Max ABSOLUTE Error: Displacement:{0:e3} , Rotation:{1:e3}",
102+
diff.Displacements.Length, diff.Rotations.Length));
103+
104+
var buf = new ValidationResult();
105+
buf.Span = span;
106+
buf.Title = "Rigid element Validation, Console beam with tetrahedron";
107+
108+
return buf;
109+
}
110+
}
111+
}
112+
}

‎BriefFiniteElementNet.Validation/GithubIssues/Issue96.cs‎

Lines changed: 45 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,45 @@
1+
using System;
2+
using System.Collections.Generic;
3+
using System.Linq;
4+
using System.Text;
5+
using System.Threading.Tasks;
6+
using BriefFiniteElementNet.Elements;
7+
8+
namespace BriefFiniteElementNet.Validation.GithubIssues
9+
{
10+
public class Issue96
11+
{
12+
public static void Run1()
13+
{
14+
var m1 = new Model();
15+
16+
m1.Nodes.Add(new Node(0, 0, 0) { Label = "n1" });
17+
18+
m1.Nodes.Add(new Node(0, 0, 6) { Label = "n2" });
19+
20+
m1.Nodes["n1"].Constraints = new Constraint(dx: DofConstraint.Fixed, dy: DofConstraint.Fixed, dz: DofConstraint.Fixed, rx: DofConstraint.Fixed, ry: DofConstraint.Fixed, rz: DofConstraint.Fixed);
21+
22+
m1.Nodes["n2"].Loads.Add(new NodalLoad(new Force(0, 0, -4000, 0, 0, 0)));
23+
24+
m1.Elements.Add(new BarElement(m1.Nodes["n1"], m1.Nodes["n2"]) { Label = "r1" });
25+
26+
(m1.Elements["r1"] as BarElement).Behavior = BarElementBehaviours.FullFrame;
27+
28+
(m1.Elements["r1"] as BarElement).Section = new BriefFiniteElementNet.Sections.UniformParametric1DSection(16.43 / 10000, 541.22 / 100000000, 44.92 / 100000000, 541.22 / 100000000 + 44.92 / 100000000);
29+
30+
(m1.Elements["r1"] as BarElement).Material = BriefFiniteElementNet.Materials.UniformIsotropicMaterial.CreateFromYoungPoisson(210 * Math.Pow(10, 9), 0.3);
31+
32+
m1.Solve();
33+
34+
var r11 = m1.Nodes[0].GetSupportReaction();
35+
36+
var r12 = m1.Nodes[1].GetSupportReaction();
37+
38+
Console.WriteLine("n0 reaction: {0}", r11);
39+
40+
Console.WriteLine("n1 reaction: {0}", r12);
41+
42+
Console.ReadKey();
43+
}
44+
}
45+
}

‎BriefFiniteElementNet/ElementHelpers/TetrahedronHelper.cs‎

Lines changed: 16 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -13,7 +13,7 @@ public class TetrahedronHelper : IElementHelper
1313
public Element TargetElement { get; set; }
1414

1515
/// <inheritdoc/>
16-
public Matrix GetBMatrixAt_new(Element targetElement, params double[] isoCoords)
16+
public Matrix GetBMatrixAt(Element targetElement, params double[] isoCoords)
1717
{
1818
//port from D3_TETRAH.m from fem_toolbox
1919

@@ -34,13 +34,15 @@ public Matrix GetBMatrixAt_new(Element targetElement, params double[] isoCoords)
3434
}
3535
}
3636

37+
//J.TransposeInPlace();
38+
3739
var v = 1 / 6.0 * J.Determinant();
3840

3941
var buf = targetElement.AllocateFromPool(6, 12);
4042

4143
{
4244
var Q = J.Inverse();//* J.Determinant();
43-
Q.Scale(J.Determinant());
45+
//Q.Scale(J.Determinant());
4446

4547
double a1 = Q[0, 1]; double b1 = Q[0, 2]; double c1 = Q[0, 3]; //double d1 = Q[0, 3];
4648
double a2 = Q[1, 1]; double b2 = Q[1, 2]; double c2 = Q[1, 3]; //double d2 = Q[1, 3];
@@ -60,7 +62,7 @@ public Matrix GetBMatrixAt_new(Element targetElement, params double[] isoCoords)
6062
return buf.AsMatrix();
6163
}
6264

63-
public Matrix GetBMatrixAt(Element targetElement, params double[] isoCoords)
65+
public Matrix GetBMatrixAt_old(Element targetElement, params double[] isoCoords)
6466
{
6567
//http://what-when-how.com/the-finite-element-method/fem-for-3d-solids-finite-element-method-part-1/
6668
//https://academic.csuohio.edu/duffy_s/CVE_512_12.pdf
@@ -201,6 +203,17 @@ public Matrix GetJMatrixAt(Element targetElement, params double[] isoCoords)
201203
/// <inheritdoc/>
202204
public Matrix CalcLocalStiffnessMatrix(Element targetElement)
203205
{
206+
//edit
207+
var b = GetBMatrixAt(targetElement,0, 0, 0);
208+
var d = GetDMatrixAt(targetElement, 0, 0, 0);
209+
210+
var k = b.Transpose() * d * b;
211+
212+
var v = GetJMatrixAt(targetElement, 0, 0, 0).Determinant() / 6.0;
213+
214+
k.Scale(v);
215+
216+
return k;
204217
return ElementHelperExtensions.CalcLocalKMatrix_Tetrahedron(this, targetElement);
205218
}
206219

0 commit comments

Comments
 (0)