BriefFiniteElementNet / BriefFiniteElement.Net

BriefFiniteElementDotNET (BFE.NET) is a library for linear-static Finite Element Method (FEM) analysis of solids and structures in .NET
GNU Lesser General Public License v3.0
154 stars 57 forks source link

ExactInternalForces incorrect with partial BarElement end releases #48

Closed NickN-Eng closed 4 years ago

NickN-Eng commented 4 years ago

Hi epsi1on,

First of all, a big thank-you for the work you have done on this project!

I was wondering if the BarElement Start/EndReleaseConditions were supposed to be working? Alternatively, it may be an issue with the internal forces (I am using GetExactInternalForceAt()).

In the problem below, I have defined a simply supported beam with 1kN/m uniformly distributed load. However, instead of fixing the Start/End NODES, I am setting the EndReleaseConditions. I have fixed the Rz of one of the beam ends so that the bar is fully constrained.

The internal forces are incorrect. The support reactions are correct.

image

I want to set the EndReleaseCondition on the bar, instead of the node constraint because I want to torsion ally restrain the beam in it's own axis. For example, if I had a diagonal beam in 3D, I cannot restrain the beam in pure torsion (local Rz) through nodal restraints - nodals restraints are in global coordinates.

Code is below:

        var l = 10;
        var model = new Model();
        model.Nodes.Add(new Node(0, 0, 0) { Label = "n0" });
        model.Nodes.Add(new Node(l, 0, 0) { Label = "n1" });
        //Fixed nodes
        model.Nodes["n0"].Constraints = Constraints.Fixed;
        model.Nodes["n1"].Constraints = Constraints.Fixed;

        var bar = new BarElement(model.Nodes["n0"], model.Nodes["n1"]) { Label = "e0" };
        //Pinned bar releases
        bar.StartReleaseCondition = Constraints.MovementFixed;
        bar.EndReleaseCondition = new Constraint(dx: DofConstraint.Released, dy: DofConstraint.Fixed, dz: DofConstraint.Fixed, rx: DofConstraint.Fixed, ry: DofConstraint.Released, rz: DofConstraint.Released);
        model.Elements.Add(bar);

        var sec = new UniformGeometric1DSection(SectionGenerator.GetISetion(0.24, 0.67, 0.01, 0.006));
        var mat = UniformIsotropicMaterial.CreateFromYoungPoisson(210e9, 0.3);

        (model.Elements["e0"] as BarElement).Material = mat;
        (model.Elements["e0"] as BarElement).Section = sec;

        var u1 = new UniformLoad(LoadCase.DefaultLoadCase, -Vector.K, 1, CoordinationSystem.Global);

        model.Elements["e0"].Loads.Add(u1);

        model.Solve_MPC();

        var n0Force = model.Nodes["n0"].GetSupportReaction();
        var n1Force = model.Nodes["n1"].GetSupportReaction();
        Console.WriteLine("support reaction of n0: {0}, n1: {1}", n0Force, n1Force);

        var elm = model.Elements[0] as BarElement;
        var fnc = new Func<double, double>(x =>
        {
            try
            {
                var xi = elm.LocalCoordsToIsoCoords(x);
                var frc = elm.GetExactInternalForceAt(xi[0]);
                return frc.My;
            }
            catch
            {

                return 0;
            }

        });

        FunctionVisualizer.VisualizeInNewWindow(fnc, 1E-6, l - 1E-6, 50);
epsi1on commented 4 years ago

I've just commented line 1028 through 1043 of file EulerBernoulliBeamHelper.cs and this issue is fixed. Not sure, but i think did add those lines to fix another issue. Anyways this is the result for your code: image

NickN-Eng commented 4 years ago

Thanks for the quick response! Have had a look at the version history and you made this update in Commit 872e35b2 07/11/2019 13:37:47 "update". I dont fully understand how GetLoadInternalForceAt() works, but as part of this update, you removed this: image ...and added: image Do you remember why you made the change? I hope commenting out lines 1028->1043 doesnt break something else!!

epsi1on commented 4 years ago

I've added the unit test for this issue:

BriefFiniteElementNet.Tests.BarElementUniformLoadsTests.LoadInternalForce_uniformload_eulerbernoullybeam_endrelease()

This way if any other changes break this issue, we'll understand.

NickN-Eng commented 4 years ago

Cheers epsi1on!