Simple finite elements matrix structural analysis using GPU. Actually supports 2D planar membrane/shell.
C++ STANDARD 14 and CUDA. No other libraries are needed.
To fully use this library, you must include the following files:
#include "fnelem/math/fematrix.cpp"
#include "fnelem/math/fematrix_utils.cpp"
#include "fnelem/math/matrix_inversion_cpu.cpp"
#include "fnelem/math/matrix_inversion_cuda.cu"
#include "fnelem/analysis/static_analysis.cpp"
#include "fnelem/model/base/model.cpp"
#include "fnelem/model/base/model_component.cpp"
#include "fnelem/model/elements/element.cpp"
#include "fnelem/model/elements/membrane.cpp"
#include "fnelem/model/loads/load.cpp"
#include "fnelem/model/loads/load_membrane_distributed.cpp"
#include "fnelem/model/loads/load_node.cpp"
#include "fnelem/model/loads/load_pattern.cpp"
#include "fnelem/model/loads/load_pattern_constant.cpp"
#include "fnelem/model/nodes/node.cpp"
#include "fnelem/model/restraints/restraint.cpp"
#include "fnelem/model/restraints/restraint_node.cpp"The model class is the base object for defining a problem. It needs the number of degrees of freedom, the dimension (actually, only 2D is supported), the nodes, elements, restraints, loads, and load patterns.
Model *model = new Model(dimension, number_degrees_of_freedom);
// Define nodes
model->add_nodes(nodes);
model->add_elements(elements);
model->add_restraints(restraints);
model->add_load_patterns(load_pattern);Nodes, elements, restraints, and load_pattern are std::vector objects that contain the library Objects.
The StaticAnalysis class performs the basic static analysis method. Using the model definition, it creates the structure stiffness matrix, and then using K*u = F, it solves the displacement vector of the nodes.
StaticAnalysis *analysis = new StaticAnalysis(model);
analysis->analyze(use_gpu);If use_gpu is true then GPU performs the main matrix inversion. After the analysis is complete, the model can save the results into a file. That file contains the main structure elements: nodes, shells, reactions, and internal forces of the elements.
model->save_results("file.txt");Example of the output:
FNELEM-GPU - Finite element structural analysis using CUDA and GPU.
v0.2.0 (21/12/2018) @ ppizarror
--------------------------------------------------------------------
Input model properties:
--------------------------------------------------------------------
Nodes:
Node count: 6
Node N1: 0.000000 0.000000
...
Elements:
Element count: 2
Membrane MEM1:
Width (2b): 2
Height (2h): 2
Thickness: 1
Elastic mod: 2000
Poisson mod: 0.2
Element nodes: N1, N3, N4, N2
...
--------------------------------------------------------------------
Analysis results:
--------------------------------------------------------------------
Node displacements:
Node N1: 0.000000 0.000000
...
Node reactions:
Node N1: 369.158879 500.000000
...
Element stresses:
Membrane MEM1:
Node N1 (-b, -h): 369.159, 500
Node N3 (+b, -h): 130.841, 0
Node N4 (+b, +h): -500, -500
Node N2 (-b, +h): 0, 0
Stress MEM1 [GLX GLY X Y SIGMAX SIGMAY SIGMAXY DISPLX DISPLY]
0 0 -1 -1 -44.3925 -221.963 -193.9250 0
...
The node class must be defined using the global coordinates. Nodes are used by the elements.
std::vector<Node *> *nodes = new std::vector<Node *>();
nodes->push_back(new Node(node_tag, coordinate_x, coorinate_y));
nodes->push_back(new Node(node_tag, coordinate_x, coordinate_y, coordinate_z));The membrane shell must be defined using four nodes, the elastic and Poisson modulus, and the thickness.
std::vector<Element *> *elements = new std::vector<Element *>();
elements->push_back(new Membrane(membrane_tag, node1, node2, node3, node4, E, poisson, thickness));The restraint class must be used when a node cannot move in any direction.
std::vector<Restraint *> *restraints = new std::vector<Restraint *>();
RestraintNode *r = new RestraintNode(restraint_tag, node);
r1->add_dofid(index);
r1->add_all(); // This method adds all degrees of freedom
restraints->push_back(r);NodeLoad is a load applied to a node; it needs the node reference and a load vector. In FNELEM-GPU matrices and vectors are defined using FEMatrix class. That class provides the following main methods:
- Add
- Subtract
- Multiply
- Inverse (Using CPU and GPU)
- Determinant
- max(), min()
- clone()
To define a vector:
FEMatrix *vector = FEMatrix_vector(size);
vector->set(index, value);
vector->get(index); To define a matrix:
FEMatrix *matrix = FEMatrix_vector(n, m);
matrix->set(i, j, value);
matrix->get(i, j); Then, node load is defined by:
// Create load vector
std::vector<Load *> *loads = new std::vector<Load *>();
// Create LoadNode
FEMatrix *load = FEMatrix_vector(size);
load->set(index, value);
...
loads->push_back(new LoadNode(load_tag, node, load));This class can apply a distributed load. It needs the membrane, the node position (1, 2, 3, 4)
std::vector<Load *> *loads = new std::vector<Load *>();
loads->push_back(new LoadMembraneDistributed(load_tag, membrane, position_initial,
position_final, load_at_initial, distance_from_initial, load_at_final, distance_from_initial));The load pattern must be used to apply the loads to the model. Only the constant load pattern is defined.
std::vector<LoadPattern *> *loadpattern = new std::vector<LoadPattern *>();
loadpattern->push_back(new LoadPatternConstant(load_tag, loads));The following code defines a simple bridge composed of N membranes under a distributed load of 100kN/m. The membranes are the same, 100cm height, 100cm width, 15cm thickness, elastic modulus of 300000kN/cm^2, Poisson modulus of 0.2. This example can also be found in test/analysis/test_static_analysis.h.
int N = 3; // Number of membranes
int b = 100; // Pier width
int h = 100; // Pier height
// N = 2
// 4 ------ 5 ------ 6
// | | |
// h | (1) | (2) | ......
// | | |
// 1 ------ 2 ------ 3 .....
// ^ b ^
//
// N = 3
// 4 ------ 5 ------ 6 ------- 8
// | | | |
// h | (1) | (2) | (3) |
// | | | |
// 1 ------ 2 ------ 3 ------- 4
// ^ b ^
// =====================================
double t = 15; // Thickness (cm)
double E = 300000; // Elastic modulus
double nu = 0.15; // Poisson modulus
// Number degrees of freedom
int gdl = N * 4;
// Create model
Model *model = new Model(2, gdl);
// Create nodes
std::vector<Node *> *nodes = new std::vector<Node *>();
int j;
for (int i = 1; i <= N + 1; i++) {
nodes->push_back(new Node("N" + std::to_string(i), b * (i - 1), 0));
}
for (int i = 1; i <= N + 1; i++) {
j = N + 1 + i;
nodes->push_back(new Node("N" + std::to_string(j), b * (i - 1), h));
}
// Add nodes to model
model->add_nodes(nodes);
// Create elements
// n4 ------------ n3
// | |
// | (i) |
// | |
// n1 ------------ n2
unsigned long n1, n2, n3, n4;
std::vector<Element *> *elements = new std::vector<Element *>();
for (unsigned long i = 0; i < N; i++) {
n1 = i;
n2 = i + 1;
n3 = N + i + 2;
n4 = N + i + 1;
elements->push_back(new Membrane("MEM" + std::to_string(i), nodes->at(n1),
nodes->at(n2), nodes->at(n3), nodes->at(n4), E, nu, t));
}
model->add_elements(elements);
// Create restraints
std::vector<Restraint *> *restraints = new std::vector<Restraint *>();
RestraintNode *r1 = new RestraintNode("R1", nodes->at(0));
RestraintNode *r2 = new RestraintNode("R2", nodes->at(static_cast<unsigned long>(N)));
r1->add_all();
r2->add_all();
restraints->push_back(r1);
restraints->push_back(r2);
model->add_restraints(restraints);
// Add distribuited load
std::vector<Load *> *loads = new std::vector<Load *>();
for (int i = 0; i < elements->size(); i++) {
Membrane *mem = dynamic_cast<Membrane *>(elements->at
(static_cast<unsigned long>(i)));
loads->push_back( new LoadMembraneDistributed("DV100kN V @" +
std::to_string(i + 1), mem, 4, 3, -100, 0, -100, 1));
}
std::vector<LoadPattern *> *loadpattern = new std::vector<LoadPattern *>();
loadpattern->push_back(new LoadPatternConstant("LOADCONSTANT", loads));
model->add_load_patterns(loadpattern);
// Create analysis
StaticAnalysis *analysis = new StaticAnalysis(model);
analysis->analyze(false);
// Save results to file
model->save_results("bridge.txt");
// Delete data
analysis->clear();
delete elements;
delete loadpattern;
delete loads;
delete restraints;
delete nodes;
delete model;
delete analysis;The results are the following:
If N number of membranes is greater then GPU analysis will be faster.
This project is licensed under MIT [https://opensource.org/licenses/MIT/]
Pablo Pizarro R. | 2018 - 2019