11,927,167 members (50,066 online)
alternative version

37.2K views
56 bookmarked
Posted

# Finite Element Method Programming in C#

, 6 Apr 2015 GPL3
 Rate this:
How to use BriefFiniteElement.NET to analyse solids and structures

## Introduction

You can check the library features from Features of BriefFiniteElement.NET

As every several days the source code is changing (because of adding new features or fixing bugs) so I suggest you to download the source code and simply build them to be sure that you have latest binaries of BriefFiniteElement.NET.

# Examples

Before starting the examples, you have to add a reference to `BriefFiniteElement.NET` project or dll file.

## Example 1: A Simple Truss

In this example, I want to analyse a simple truss with 4 members as shown in picture.

All members sections are the same, a square steel section with dimension of 3 cm. So the properties of members will be:

E = 210 GPa = 210*10^9 Pa = 210*10^9 N/M^2

A = 0.03m * 0.03m = 9 * 10^-4 m^2

We should do 5 steps before we solve the model:

1. Create Model, Members and Nodes.
2. Add the Nodes and Elements to Model
3. Assign geometrical and mechanical properties to Elements
4. Assign Constraints to Nodes (fix the DoF s)

And finally solve model with `Model.Solve() `method and then extract analysis results like support reactions or member internal forces or nodal deflections.

### Creating Model, Members and Nodes

#### Creating Model

We should create a Finite Element model first and then add members and nodes to it:

```// Initiating Model, Nodes and Members
var model = new Model();```

#### Creating Nodes

We should create nodes like this. In `BriefFiniteElement.NET`, every node and element have a property of type `string `named `Label `and another one named `Tag `both of which are inherited from `BriefFiniteElementNet.StructurePart`. In every `Model`, `Label `of every member should be unique among all members (both Nodes and Elements) unless the `Label `be equal to `null `which is by default. In the below code, we are creating 5 nodes of truss and assigning a unique Label to each one.

```var n1 = new Node(1, 1, 0);
n1.Label = "n1";//Set a unique label for node
var n2 = new Node(-1, 1, 0) {Label = "n2"};//using object initializer for assigning Label
var n3 = new Node(1, -1, 0) {Label = "n3"};
var n4 = new Node(-1, -1, 0) {Label = "n4"};
var n5 = new Node(0, 0, 1) {Label = "n5"};```

#### Creating Elements

Then we have to create the elements. In `BriefFiniteElement.NET`, the `TrussElement `class represents a truss element in 3D.

```var e1 = new TrussElement(n1, n5) {Label = "e1"};
var e2 = new TrussElement(n2, n5) {Label = "e2"};
var e3 = new TrussElement(n3, n5) {Label = "e3"};
var e4 = new TrussElement(n4, n5) {Label = "e4"};```

### 2- Adding Nodes and Elements to Model

You can simply add the elements and nodes we created into the `Model`. `Model `has two members Elements and Nodes which both represents an `IList<T> `of nodes and members, plus an `AddRange `method.

```model.Nodes.AddRange(n1, n2, n3, n4, n5);

Please note that if `Node `or `Element`’s `Label `property is something else than `null`, then it should be unique among all nodes and elements.

### 3- Assigning geometrical and mechanical properties to Elements

As elastic module for all members equals to 210 GPa and area of all members equals to 0.0009 m^2 we can set the element properties like this:

```e1.A = e2.A = e3.A = e4.A = 9e-4;
e1.E = e2.E = e3.E = e4.E = 210e9;```

### 4- Assigning Constraint to Nodes to fix the DoFs

Now, we should make some DoFs of structure fix in order to make analysis logically possible.

In `BriefFiniteElement.NET`, every node has 6 degree of freedom: X, Y, and Z rotations and X, Y, and Z translations. For a every truss model, we have to fix rotational DoFs for each Node (X,Y and Z rotation). Also the nodes 1 to 4 are also movement fixed, then nodes 1 to 4 should be totally fixed and node 5 should be rotation fixed. In `BriefFiniteElement.NET`, a `struct `named `Constraint `represents a constraint that is applicable to a 6 DoF node, it have Dx, Dy, Dz, Rx, Ry and Rz properties of type `DofConstraint `which is an `enum `and have two possible values 0 (Released) and 1 (Fixed). For making work easier, the `Constraint struct `has some predefined `Constraint`s in its `static `properties for example `Constraint.Fixed `or `Constraint.Free`. Here is more detailed information:

 Property name Description `Fixed` All 6 DoFs are fixed `Free` All 6 DoFs are released `MovementFixed` 3 translation DoFs are fixed and 3 rotation DoFs are released `RotationFixed` 3 translation DoFs are released and 3 rotation DoFs are fixed

We can fix DoFs of nodes 1 to 4 like this:

`n1.Constraints = n2.Constraints = n3.Constraints = n4.Constraints = new Constraint(dx:DofConstraint.Fixed, dy:DofConstraint.Fixed, dz:DofConstraint.Fixed, rx:DofConstraint.Fixed, ry:DofConstraint.Fixed, rz:DofConstraint.Fixed);`

or:

`n1.Constraints = n2.Constraints = n3.Constraints = n4.Constraints = Constraint.Fixed`

and should fix the rotational DoFs of node 5:

`n5.Constraints = Constraint.RotationFixed`

### 5- Assigning Loads to Nodes

In `BriefFiniteElement.NET`, there is a `struct `named `Force `which represent a concentrated force in 3D space which contains of 3 force components in X, Y and Z directions and three moment components in X, Y and Z directions. It have 6 double properties named Fx, Fy, Fz, Mx, My and Mz that are representing the load components. There are also two properties of type `Vector `for this `struct `named `Forces `and `Moments`. On setting or getting, they will use the Fx, Fy, Fz, Mx, My and Mz to perform operations:

```/// <summary>
/// Gets or sets the forces.
/// </summary>
/// <value>
/// The forces as a <see cref="Vector"/>.
/// </value>
public Vector Forces
{
get
{
return new Vector(fx,fy,fz);
}

set
{
this.fx = value.X;
this.fy = value.Y;
this.fz = value.Z;
}
}```

Same is with `Moments `property. The `Forces `and `Moments `property do not actually store values in something other than 6 obvious properties.

As `LoadCase `and `LoadCombination `concepts are supported in `BriefFiniteElement.NET`, every `Load `should have a `LoadCase`. A `LoadCase `is simply a `struct `that has two properties: `CaseName `with `string `type and `LoadType `with `LoadType `type which is an `enum `and has some possible values:

```public enum LoadType
{
Default = 0,
Live,
Snow,
Wind,
Quake,
Crane,
Other
}```

The `LoadType.Default `is a load type that is created for built in usage in library and it do not meant to have meaning like `Dead`, `Live`, etc. The `LoadCase struct `has a `static `property named `LoadCase.DefaultLoadCase`:

```/// <summary>
/// Gets the default load case.
/// </summary>
/// <value>
/// </value>
/// <remarks>
{
get { return new LoadCase(); }
} ```

Which represents a `LoadCase `with `LoadType `of `Default `and `CaseName `of `null`. We will call such a `LoadCase `as `DefaultLoadCase`. For simplicity of usage in `BriefFiniteElement.NET `everywhere that you’ll prompt for a `LoadCase`, if you do not provide a `LoadCase `then the `LoadCase `is assumed `DefualtLoadCase `by the library. For example, when you want to assign a load to a node, you should provide a `LoadCase `for it, like this:

`var load = new NodalLoad(new Force(0, 0, -1000, 0, 0, 0), new LoadCase("Case1",LoadType.Dead));`

but if you do not provide the `LoadCase `in the above code like this:

`var load = new NodalLoad(new Force(0, 0, -1000, 0, 0, 0));`

then the load case will be assumed `DefaultLoadCase `by the library.

Ok, next we have to add 1KN load to node 5 like this, will do it with `DefaultLoadCase`:

```var force = new Force(0, 0, -1000, 0, 0, 0);

And finally solve the model with `model.Solve() `method. Actually solving the model is done in two stages:

1. First stage is creating stiffness matrix and factorizing stiffness matrix which will take majority of time for analyzing
2. Second phase is analyzing structure against each load case which takes much less time against first stage (say for example 13 sec for first stage and 0.5 sec for second stage).

First stage is done in `model.Solve() `method and second stage will done if they’ll be need to.

There are loads with different `LoadCase`s that are applied to the Nodes and Elements. So the `Node.GetSupportReaction() `method have an overload which gets a `LoadCombination `and returns the support reactions based on the load combination. `LoadCombination `has a `static `property named `LoadCombination.DefaultLoadCombination `which has only one `LoadCase `in it (the `DefaultLoadCase`) with factor of 1.0. also everywhere that you should provide a `LoadCombination`, if you do not provide any, then `DefaultLoadCombination `will be considered by library. I’ve used `DefaultLoadCase `and `DefaultLoadCombination `in library to make working with library easier for people who are not familiar with load case and load combination stuff.

For getting the support reaction for the truss, we can simply call `Node.GetSupportReaction() `to get support reaction for every node:

```Force r1 = n1.GetSupportReaction();
Force r2 = n2.GetSupportReaction();
Force r3 = n3.GetSupportReaction();
Force r4 = n4.GetSupportReaction();```

The plus operator is overloaded for `Force struct`, so we can check the sum of support reactions:

`var rt = r1 + r2 + r3 + r4;//shows the Fz=1000 and Fx=Fy=Mx=My=Mz=0.0`

The forces (Fx, Fy and Fz) amount should be equal to sum of external loads and direction should be opposite to external loads to satisfy the structure static equilibrium equations.

### All Codes Together

This is all codes above for truss example.

```private static void Example1()
{
Console.WriteLine("Example 1: Simple 3D truss with four members");

// Initiating Model, Nodes and Members
var model = new Model();

var n1 = new Node(1, 1, 0);
n1.Label = "n1";//Set a unique label for node
var n2 = new Node(-1, 1, 0) {Label = "n2"};//using object initializer for assigning Label
var n3 = new Node(1, -1, 0) {Label = "n3"};
var n4 = new Node(-1, -1, 0) {Label = "n4"};
var n5 = new Node(0, 0, 1) {Label = "n5"};

var e1 = new TrussElement(n1, n5) {Label = "e1"};
var e2 = new TrussElement(n2, n5) {Label = "e2"};
var e3 = new TrussElement(n3, n5) {Label = "e3"};
var e4 = new TrussElement(n4, n5) {Label = "e4"};
//Note: labels for all members should be unique, else you will receive InvalidLabelException when adding it to model

e1.A = e2.A = e3.A = e4.A = 9e-4;
e1.E = e2.E = e3.E = e4.E = 210e9;

//Applying restrains

n1.Constraints = n2.Constraints = n3.Constraints = n4.Constraints = Constraint.Fixed;
n5.Constraints = Constraint.RotationFixed;

var force = new Force(0, 0, -1000, 0, 0, 0);

model.Solve();

var r1 = n1.GetSupportReaction();
var r2 = n2.GetSupportReaction();
var r3 = n3.GetSupportReaction();
var r4 = n4.GetSupportReaction();

var rt = r1 + r2 + r3 + r4;//shows the Fz=1000 and Fx=Fy=Mx=My=Mz=0.0
}
```

In this example I’m going to analyze a single sloped steel frame which its elements are under uniform distributed loads as shown in image.

The purpose is to find the internal forces of members based on shown loads.

### Creating members and nodes

First thing we are going to do is to create the members and assign the geometrical and mechanical properties to them

*Note: in next code, two columns are ‘e1’ and ‘e4’ and two beams are ‘e2’ and ‘e3’.

```var model = new Model();

var n1 = new Node(-10, 0, 0);
var n2 = new Node(-10, 0, 6);
var n3 = new Node(0, 0, 8);
var n4 = new Node(10, 0, 6);
var n5 = new Node(10, 0, 0);

var secAA = SectionGenerator.GetISetion(0.24, 0.67, 0.01, 0.006);
var secBB = SectionGenerator.GetISetion(0.24, 0.52, 0.01, 0.006);

var e1 = new FrameElement2Node(n1, n2);
var e2 = new FrameElement2Node(n2, n3);
var e3 = new FrameElement2Node(n3, n4);
var e4 = new FrameElement2Node(n4, n5);

e1.Geometry = e4.Geometry = secAA;
e2.Geometry = e3.Geometry = secBB;

e1.E = e2.E = e3.E = e4.E = 210e9;

e1.UseOverridedProperties = e2.UseOverridedProperties = e3.UseOverridedProperties = e4.UseOverridedProperties = false;```

In briefFiniteElement.NET, there is two ways to define section for a frame member. First is to define each of geometrical parameters (like area and second moments of area) by settings the FrameElement2Node.A or FrameElement2Node.Iy or FrameElement2Node.Iz, another way to do this is to define a polygon as section for frame element and library will automatically calculate the gemetrical parameters of section like area and second moments of area. By default library assumes that you are defining geometrical parameters by settings FrameElement2Node.A, Iy or Iz. if you want to provide a polygon as section, you should set the FrameElement2Node.Geometry property (which its type is PolygonYz) and right after that you should set the FrameElement2Node.UseOverridedProperties to false in order to say the library to use Geometry instead of A, Iy, Iz etc. There is also a useful static SectionGenerator class which provides some predefined sections like I sections or rectangular sections.
After creating elements and nodes, we can visualize the model in debug mode like image below.

*Note: elements with UseOverridedProperties = true are shown with square sections (dimentsion of section automatically tunes for better visualization of elements) but elements with UseOverridedProperties = false will be shown whith their real section with real dimestion.

Next we have to set boundary conditions or supports. As this problem is 2D (in x-z plane) we have to fix all movement in Y direction and also fix all rotations in x and z directions. As n1 and n5 are hinged supports we also have to fix the movement of n1 and n5 in both x and z directions too.

```n1.Constraints = n2.Constraints = n3.Constraints = n4.Constraints = n5.Constraints = Constraint.FixedDy & Constraint.FixedRx & Constraint.FixedRz;//Dy fixed and Rx fixed and Rz fixed

n1.Constraints = n1.Constraints & Constraint.MovementFixed;
n5.Constraints = n5.Constraints & Constraint.MovementFixed;```

Just have to note that the ‘&’ operator between two Constraint objects work like this: if both or one of Dx s are fixed then resut Dx is fixed, else Dx is released. It’s like Boolean operation while 'fixed' means '1' and 'free' means '0' !

There are two distributed loads which are applying to model. Both are 1 Tonf/m which is equals to 10000 N/m but there is a difference as is illustrated before. The load which is applying to e2 element is in vertical direction but the load is applying to e3 element is perpendicular to e3 element. So it can be said the load on element e2 is a uniform load in GLOBAL z direction with amount of -10000 N/m and load which is applying to e3 is in element LOCAL z direction of coordination system and amount of -10000 N/m.

```var ll = new UniformLoad1D(-10000, LoadDirection.Z, CoordinationSystem.Global);

### Checking for errors in model

In this example if we do not set e1.Geometry property, will get an error when equation system is solving because of structure is not stable. in this case finding the error maybe seems simple but its not always like this. for this reason there is a Trace property for Model class which will let user to see status of model while its working. For example on calling Model.Solve() for first time, several separated task are being done. first is calculating stiffness matrix of each element, next they should be assembeled into global stiffness matrix KG, then kff should extract from KG and finally equation system that contains kff should be solved. User do not see the information but information about these datas are possible to achieve from Model.Trace. It like System.Diagnostics.Trace functionality. The Trace member have TraceListeners property which is a list of ITraceListener objects. once an ITraceListener object is added to this list, it will receive trace info of model. There is a ready ITraceListener class named WpfTraceListener in BriefFiniteElementNet.Controls project that can be used for seeing the information. here is example:

```var wnd = WpfTraceListener.CreateModelTrace(model);
model.CheckForErrors();
wnd.ShowDialog();```

the first line creates a WpfTraceListener object and add to Model.Trace.TraceListeners list. model.CheckForErrors () will check model for usual (obvious) errors and if any found, will be written into Model.Trace so we can see it our window. and finally wnd.ShowDialog() will shows the window. if for example we do not set Geometry property for e2 element, the will see this:

As you can see there is a table and a record inside of it. First field is TimeStamp which is showing the time, next is a brief message about the record, next is Error ID (in case that it be error or warning) next is Level. Level have 3 possible values:

* Error

* Warning

* Info

For example of Level, look at next pictures. Next is TargetField whici is showing the target of record, in our case it is e2 element and it have a link for more information (in this case it redirect to error list page) ad you can see description about error number MA10010 in the link.

If model do not contains error and we call wnd.ShowDialog() after call to model.Solve() like this:

```var wnd = WpfTraceListener.CreateModelTrace(model);
model.CheckForErrors();
model.Solve();
wnd.ShowDialog();```

Will see something like this:

### Analyzing the model and getting internal forces

For analyzing the model we should simply call the Model.Solve() Method. After model is solved we can get the internal force of members in any position of member's length with calling GetInternalForce method of FrameElement2Node class.
Following is the moment diagram based on output of FrameElement2Node.GetInternalForce() for this model.

## Points of Interest

At this time, the project have some features that are listed on Features of BriefFiniteElement.NET. There are also a set of features I planned to include in this library which I am working.

## Support

If you have questions about the library or you've found a bug or have any idea/feature request, please let me know via DISCUSSION or ISSUEs on Codeplex, or mail me.

## Contribute

If you have good knowledge in linear fem, look at HERE and help us solve our problems.

## TODO

• Example 1: Simple Truss (Done)
• Example 2: Element Distributed Load (Done)
• Example 4: Analyzing grid and checking result for statically equilibrium.

## History

• 9 July 2014 - First version

## Share

 Software Developer N/A Iran (Islamic Republic Of)
I'm a civil engineer with several years of experiences in C# and WPF.

## You may also be interested in...

 First PrevNext
 My vote of 5 Philip Liebscher15-Apr-15 7:23 Philip Liebscher 15-Apr-15 7:23
 Linear Load Member 113597687-Mar-15 13:34 Member 11359768 7-Mar-15 13:34
 Re: Linear Load Ehsan.MA8-Mar-15 5:41 Ehsan.MA 8-Mar-15 5:41
 Re: Linear Load Member 113597688-Mar-15 9:45 Member 11359768 8-Mar-15 9:45
 GetGlobalEquivalentNodalLoads Member 1133031523-Dec-14 5:33 Member 11330315 23-Dec-14 5:33
 Re: GetGlobalEquivalentNodalLoads Ehsan.MA23-Dec-14 5:51 Ehsan.MA 23-Dec-14 5:51
 Re: GetGlobalEquivalentNodalLoads Member 1133031524-Dec-14 1:45 Member 11330315 24-Dec-14 1:45
 Elastic line Member 1133031522-Dec-14 23:31 Member 11330315 22-Dec-14 23:31
 Re: Elastic line Ehsan.MA22-Dec-14 23:48 Ehsan.MA 22-Dec-14 23:48
 Re: Elastic line Member 1133031523-Dec-14 0:12 Member 11330315 23-Dec-14 0:12
 Re: Elastic line Ehsan.MA23-Dec-14 0:15 Ehsan.MA 23-Dec-14 0:15
 About the section property "j" Member 1119568626-Nov-14 3:21 Member 11195686 26-Nov-14 3:21
 Re: About the section property "j" Ehsan.MA26-Nov-14 3:49 Ehsan.MA 26-Nov-14 3:49
 Re: About the section property "j" Member 1119568626-Nov-14 19:10 Member 11195686 26-Nov-14 19:10
 Re: About the section property "j" Ehsan.MA26-Nov-14 19:47 Ehsan.MA 26-Nov-14 19:47
 Re: About the section property "j" Ehsan.MA26-Nov-14 4:23 Ehsan.MA 26-Nov-14 4:23
 Meaning of bitwise & operator Hiresdev3-Nov-14 19:28 Hiresdev 3-Nov-14 19:28
 Re: Meaning of bitwise & operator Ehsan.MA3-Nov-14 21:38 Ehsan.MA 3-Nov-14 21:38
 Re: Meaning of bitwise & operator Hiresdev5-Nov-14 1:47 Hiresdev 5-Nov-14 1:47
 Re: Meaning of bitwise & operator Ehsan.MA5-Nov-14 1:51 Ehsan.MA 5-Nov-14 1:51
 Re: Meaning of bitwise & operator Hiresdev5-Nov-14 2:46 Hiresdev 5-Nov-14 2:46
 Re: Meaning of bitwise & operator Ehsan.MA6-Nov-14 7:03 Ehsan.MA 6-Nov-14 7:03
 Re: Meaning of bitwise & operator Hiresdev6-Nov-14 18:58 Hiresdev 6-Nov-14 18:58
 Re: Meaning of bitwise & operator Hiresdev7-Nov-14 2:48 Hiresdev 7-Nov-14 2:48
 Re: Meaning of bitwise & operator Ehsan.MA7-Nov-14 5:40 Ehsan.MA 7-Nov-14 5:40
 Boundary Conditions for Frame Example Hiresdev2-Nov-14 3:21 Hiresdev 2-Nov-14 3:21
 Re: Boundary Conditions for Frame Example Ehsan.MA2-Nov-14 3:34 Ehsan.MA 2-Nov-14 3:34
 I have a problem with constraints... Deep Cloud21-Sep-14 10:33 Deep Cloud 21-Sep-14 10:33
 Re: I have a problem with constraints... Ehsan.MA21-Sep-14 10:54 Ehsan.MA 21-Sep-14 10:54
 Re: I have a problem with constraints... Deep Cloud21-Sep-14 11:06 Deep Cloud 21-Sep-14 11:06
 Re: I have a problem with constraints... Ehsan.MA21-Sep-14 11:10 Ehsan.MA 21-Sep-14 11:10
 Re: I have a problem with constraints... Deep Cloud21-Sep-14 11:33 Deep Cloud 21-Sep-14 11:33
 Re: I have a problem with constraints... Ehsan.MA21-Sep-14 12:08 Ehsan.MA 21-Sep-14 12:08
 Re: I have a problem with constraints... Deep Cloud21-Sep-14 11:46 Deep Cloud 21-Sep-14 11:46
 Re: I have a problem with constraints... Ehsan.MA21-Sep-14 11:53 Ehsan.MA 21-Sep-14 11:53
 Re: I have a problem with constraints... Deep Cloud21-Sep-14 12:04 Deep Cloud 21-Sep-14 12:04
 Re: I have a problem with constraints... Deep Cloud21-Sep-14 12:08 Deep Cloud 21-Sep-14 12:08
 Re: I have a problem with constraints... Ehsan.MA21-Sep-14 12:16 Ehsan.MA 21-Sep-14 12:16
 Re: I have a problem with constraints... Deep Cloud21-Sep-14 12:37 Deep Cloud 21-Sep-14 12:37
 Re: I have a problem with constraints... Deep Cloud21-Sep-14 12:15 Deep Cloud 21-Sep-14 12:15
 Can't find model.Nodes.AddRange Member 1103117223-Aug-14 3:41 Member 11031172 23-Aug-14 3:41
 Re: Can't find model.Nodes.AddRange Ehsan.MA23-Aug-14 7:50 Ehsan.MA 23-Aug-14 7:50
 Congratulations. Member 1055260618-Aug-14 2:55 Member 10552606 18-Aug-14 2:55
 Re: Congratulations. Ehsan.MA19-Aug-14 9:43 Ehsan.MA 19-Aug-14 9:43
 My vote of 5 loPetS15-Jul-14 5:41 loPetS 15-Jul-14 5:41
 Re: My vote of 5 Ehsan.MA15-Jul-14 9:53 Ehsan.MA 15-Jul-14 9:53
 My vote of 4 CatchExAs14-Jul-14 12:35 CatchExAs 14-Jul-14 12:35
 Re: My vote of 4 Ehsan.MA15-Jul-14 0:57 Ehsan.MA 15-Jul-14 0:57