Skip to content
Snippets Groups Projects
Commit d829d8ef authored by Maltez Cavalheiro Kévin's avatar Maltez Cavalheiro Kévin
Browse files

first implementation of openmp

parent 44e77f24
No related branches found
No related tags found
1 merge request!6merge kevin branch
Expand all Hide whitespace changes
Inline Side-by-side
Showing
with 1348 additions and 10 deletions
BEM.xcodeproj/project.pbxproj
+ 78
− 0
View file @ d829d8ef
......@@ -8,6 +8,9 @@
/* Begin PBXBuildFile section */
0A591E4927F9CE1700B1D4E1 /* main.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 0A591E4827F9CE1700B1D4E1 /* main.cpp */; };
0A591E5527F9D06200B1D4E1 /* functions.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 0A591E5027F9D06200B1D4E1 /* functions.cpp */; };
0A591E5627F9D06200B1D4E1 /* code_louis_function.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 0A591E5227F9D06200B1D4E1 /* code_louis_function.cpp */; };
0A591E5727F9D06200B1D4E1 /* code_louis.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 0A591E5427F9D06200B1D4E1 /* code_louis.cpp */; };
/* End PBXBuildFile section */
/* Begin PBXCopyFilesBuildPhase section */
......@@ -25,6 +28,19 @@
/* Begin PBXFileReference section */
0A591E4527F9CE1700B1D4E1 /* BEM */ = {isa = PBXFileReference; explicitFileType = "compiled.mach-o.executable"; includeInIndex = 0; path = BEM; sourceTree = BUILT_PRODUCTS_DIR; };
0A591E4827F9CE1700B1D4E1 /* main.cpp */ = {isa = PBXFileReference; lastKnownFileType = sourcecode.cpp.cpp; path = main.cpp; sourceTree = "<group>"; };
0A591E4F27F9D06200B1D4E1 /* code_louis.hpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; path = code_louis.hpp; sourceTree = "<group>"; };
0A591E5027F9D06200B1D4E1 /* functions.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = functions.cpp; sourceTree = "<group>"; };
0A591E5127F9D06200B1D4E1 /* functions.hpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; path = functions.hpp; sourceTree = "<group>"; };
0A591E5227F9D06200B1D4E1 /* code_louis_function.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = code_louis_function.cpp; sourceTree = "<group>"; };
0A591E5327F9D06200B1D4E1 /* code_louis_function.hpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; path = code_louis_function.hpp; sourceTree = "<group>"; };
0A591E5427F9D06200B1D4E1 /* code_louis.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = code_louis.cpp; sourceTree = "<group>"; };
0A591E5827F9D26900B1D4E1 /* hole.geo */ = {isa = PBXFileReference; lastKnownFileType = text; path = hole.geo; sourceTree = "<group>"; };
0A591E5927F9D26900B1D4E1 /* complex_validation.geo */ = {isa = PBXFileReference; lastKnownFileType = text; path = complex_validation.geo; sourceTree = "<group>"; };
0A591E5A27F9D26900B1D4E1 /* my_geo.geo */ = {isa = PBXFileReference; lastKnownFileType = text; path = my_geo.geo; sourceTree = "<group>"; };
0A591E5B27F9D26900B1D4E1 /* l_shape.geo */ = {isa = PBXFileReference; lastKnownFileType = text; path = l_shape.geo; sourceTree = "<group>"; };
0A591E5C27F9D26900B1D4E1 /* rectangle.geo */ = {isa = PBXFileReference; lastKnownFileType = text; path = rectangle.geo; sourceTree = "<group>"; };
0A591E5D27F9D26900B1D4E1 /* impossible.geo */ = {isa = PBXFileReference; lastKnownFileType = text; path = impossible.geo; sourceTree = "<group>"; };
0A591E5E27F9D26900B1D4E1 /* doubleHole.geo */ = {isa = PBXFileReference; lastKnownFileType = text; path = doubleHole.geo; sourceTree = "<group>"; };
/* End PBXFileReference section */
/* Begin PBXFrameworksBuildPhase section */
......@@ -58,6 +74,19 @@
isa = PBXGroup;
children = (
0A591E4827F9CE1700B1D4E1 /* main.cpp */,
0A591E5027F9D06200B1D4E1 /* functions.cpp */,
0A591E5127F9D06200B1D4E1 /* functions.hpp */,
0A591E4F27F9D06200B1D4E1 /* code_louis.hpp */,
0A591E5427F9D06200B1D4E1 /* code_louis.cpp */,
0A591E5227F9D06200B1D4E1 /* code_louis_function.cpp */,
0A591E5327F9D06200B1D4E1 /* code_louis_function.hpp */,
0A591E5927F9D26900B1D4E1 /* complex_validation.geo */,
0A591E5E27F9D26900B1D4E1 /* doubleHole.geo */,
0A591E5827F9D26900B1D4E1 /* hole.geo */,
0A591E5D27F9D26900B1D4E1 /* impossible.geo */,
0A591E5B27F9D26900B1D4E1 /* l_shape.geo */,
0A591E5A27F9D26900B1D4E1 /* my_geo.geo */,
0A591E5C27F9D26900B1D4E1 /* rectangle.geo */,
);
path = BEM;
sourceTree = "<group>";
......@@ -119,7 +148,10 @@
isa = PBXSourcesBuildPhase;
buildActionMask = 2147483647;
files = (
0A591E5627F9D06200B1D4E1 /* code_louis_function.cpp in Sources */,
0A591E5727F9D06200B1D4E1 /* code_louis.cpp in Sources */,
0A591E4927F9CE1700B1D4E1 /* main.cpp in Sources */,
0A591E5527F9D06200B1D4E1 /* functions.cpp in Sources */,
);
runOnlyForDeploymentPostprocessing = 0;
};
......@@ -239,6 +271,29 @@
isa = XCBuildConfiguration;
buildSettings = {
CODE_SIGN_STYLE = Automatic;
HEADER_SEARCH_PATHS = (
"/usr/local/include/gmsh-4.9.5-MacOSX-sdk/include",
/usr/local/eigen/3.4.0_1/include/eigen3/,
/usr/local/opt/libomp/include,
);
LIBRARY_SEARCH_PATHS = (
"/usr/local/lib/gmsh-4.9.5-MacOSX-sdk/lib",
/usr/local/eigen/3.4.0_1/include/eigen3/,
/usr/local/opt/libomp/lib,
);
OTHER_CFLAGS = (
"-Xpreprocessor",
"-fopenmp",
);
OTHER_LDFLAGS = (
"-Xclang",
"-fopenmp",
"-O2",
"-g",
"-lgmsh",
"-DNDEBUG",
"-lomp",
);
PRODUCT_NAME = "$(TARGET_NAME)";
};
name = Debug;
......@@ -247,6 +302,29 @@
isa = XCBuildConfiguration;
buildSettings = {
CODE_SIGN_STYLE = Automatic;
HEADER_SEARCH_PATHS = (
"/usr/local/include/gmsh-4.9.5-MacOSX-sdk/include",
/usr/local/eigen/3.4.0_1/include/eigen3/,
/usr/local/opt/libomp/include,
);
LIBRARY_SEARCH_PATHS = (
"/usr/local/lib/gmsh-4.9.5-MacOSX-sdk/lib",
/usr/local/eigen/3.4.0_1/include/eigen3/,
/usr/local/opt/libomp/lib,
);
OTHER_CFLAGS = (
"-Xpreprocessor",
"-fopenmp",
);
OTHER_LDFLAGS = (
"-Xclang",
"-fopenmp",
"-O2",
"-g",
"-lgmsh",
"-DNDEBUG",
"-lomp",
);
PRODUCT_NAME = "$(TARGET_NAME)";
};
name = Release;
......
BEM.xcodeproj/xcuserdata/Kevin.xcuserdatad/xcschemes/xcschememanagement.plist
+ 8
− 0
View file @ d829d8ef
......@@ -10,5 +10,13 @@
<integer>0</integer>
</dict>
</dict>
<key>SuppressBuildableAutocreation</key>
<dict>
<key>0A591E4427F9CE1700B1D4E1</key>
<dict>
<key>primary</key>
<true/>
</dict>
</dict>
</dict>
</plist>
BEM/functions.cpp 0 → 100644
+ 927
− 0
View file @ d829d8ef
This diff is collapsed.
BEM/functions.hpp 0 → 100644
+ 45
− 0
View file @ d829d8ef
#include <gmsh.h>
#include <iostream>
#include <map>
#include <sstream>
#include <Eigen/Dense>
#include <cmath>
#include <iomanip>
#include <stdio.h>
#include <chrono>
#include <omp.h>
using namespace std;
using namespace gmsh;
using namespace Eigen;
// Each element has all its information gathered as in this structure.
struct elementStruct{
int tag; // tag of the element
string bcType; // "dirichlet" or "neumann" (at the end, should replace this by an int or boolean)
double bcValue [2]; // values of [u, u_n]
size_t nodeTags [2]; // tags of the [first, second] nodes of the element
size_t midNodeTags;
double length;
double x1, y1, x2, y2;
double midX, midY;
double deltaX, deltaY;
};
int solverBEMFun(int argc, char **argv, std::vector<size_t> &finalElementTags, std::vector<double> &electrostaticPressure);
void dispNodeCoordFun(const vector<size_t> &allNodeTags, const vector<double> &allCoord);
void dispMatricesFun(const MatrixXd &A, const MatrixXd &B, const MatrixXd &c, const MatrixXd &b, const MatrixXd &x);
void dispValuesFun(const vector<elementStruct> &elementVector);
void getBC(map<std::size_t, int> allNodeIndex, const vector<double> &allCoord, vector<elementStruct> &elementVector, int &dirichletCount, const bool dispHierarchy, vector<size_t> &NodeTags2d_domain, vector<double> &NodeCoord2d_domain);
void computeGeometricParam(elementStruct &element, map<std::size_t, int> allNodeIndex, const vector<double> &allCoord);
double computeH(const double &x, const double &y, const elementStruct &element);
double computeG(const double &x, const double &y, const elementStruct &element, const string &integrationType);
void computeData(const vector<elementStruct> &elementVector, const vector<int> &elementTypes, const vector<vector<size_t>> &elementTags, vector<size_t> &tags, vector<vector<double>> &data, const string integrationType);
void displayData(const string viewName, const string modelName, const string fileName, const vector<size_t> &tags, const vector<vector<double>> &data);
void compute_internal_phi(vector<double> &allCoord,vector<size_t> &NodeTags2d_domain, vector<elementStruct> &elementVector, map<int,double> &nodeMapDomain, const string integrationType);
void element_node_data(vector<vector<vector<double>>> &data, vector<vector<size_t>> &elementTags2d, vector<int> &elementTypes2d, vector<vector<size_t>> &nodeTags2d, vector<elementStruct> &elementVector, map<int, double> &nodeMapDomain, bool &show_element_node_data);
void disp_element_node_data(vector<int> &nodes_commun_tag, vector<int> &elements_commun_tag);
BEM/hole.geo 0 → 100644
+ 60
− 0
View file @ d829d8ef
Lx_in = 1.0;
Ly_in = 0.1;
cx = 0.0;
cy = -0.6;
Lx_out = 2.0;
Ly_out = 1.5;
phi_in = 200.;
phi_out = 100.;
nx_in = 40;
ny_in = 4;
nx_out = 80;
ny_out = 60;
Point(1) = {-Lx_out/2, -Ly_out/2, 0, 1.0};
Point(2) = {Lx_out/2, -Ly_out/2, 0, 1.0};
Point(3) = {Lx_out/2, Ly_out/2, 0, 1.0};
Point(4) = {-Lx_out/2, Ly_out/2, 0, 1.0};
Point(5) = {-Lx_in/2 + cx, -Ly_in/2 + cy, 0, 1.0};
Point(6) = {-Lx_in/2 + cx, Ly_in/2 + cy, 0, 1.0};
Point(7) = {Lx_in/2 + cx, Ly_in/2 + cy, 0, 1.0};
Point(8) = {Lx_in/2 + cx, -Ly_in/2 + cy, 0, 1.0};
Line(1) = {1, 2};
Line(2) = {2, 3};
Line(3) = {3, 4};
Line(4) = {4, 1};
Line(5) = {5, 6};
Line(6) = {6, 7};
Line(7) = {7, 8};
Line(8) = {8, 5};
Curve Loop(1) = {4, 1, 2, 3};
Curve Loop(2) = {8, 5, 6, 7};
Plane Surface(1) = {1, 2};
Transfinite Curve {1} = nx_out + 1 Using Progression 1;
Transfinite Curve {2} = ny_out + 1 Using Progression 1;
Transfinite Curve {3} = nx_out + 1 Using Progression 1;
Transfinite Curve {4} = ny_out + 1 Using Progression 1;
Transfinite Curve {5} = ny_in + 1 Using Progression 1;
Transfinite Curve {6} = nx_in + 1 Using Progression 1;
Transfinite Curve {7} = ny_in + 1 Using Progression 1;
Transfinite Curve {8} = nx_in + 1 Using Progression 1;
//Recombine Surface {1};
Physical Curve("inside", 5) = {5, 6, 7, 8};
Physical Curve("outside", 6) = {1, 2, 3, 4};
SetNumber("Boundary Conditions/inside/dirichlet", phi_in);
SetNumber("Boundary Conditions/outside/dirichlet", phi_out);
\ No newline at end of file
BEM/impossible.geo 0 → 100644
+ 35
− 0
View file @ d829d8ef
L = 1.0;
phi_bottom = 200.;
phi_sides = 100.;
n = 100;
Point(1) = {0, 0, 0, 1.0};
Point(2) = {L, 0, 0, 1.0};
Point(3) = {L, L, 0, 1.0};
Point(4) = {0, L, 0, 1.0};
Line(1) = {1, 2};
Line(2) = {2, 3};
Line(3) = {3, 4};
Line(4) = {4, 1};
Curve Loop(1) = {4, 1, 2, 3};
Plane Surface(1) = {1};
Transfinite Curve {3, 1} = n+1 Using Progression 1;
Transfinite Curve {4, 2} = n+1 Using Progression 1;
Transfinite Surface {1};
Recombine Surface {1}; // quads instead of triangles
//Mesh.ElementOrder = 2;
Physical Curve("left_edge", 5) = {4};
Physical Curve("bottom_edge", 6) = {1};
Physical Curve("right_edge", 7) = {2};
Physical Curve("top_edge", 8) = {3};
// additional parameters given to the solver
SetNumber("Boundary Conditions/left_edge/dirichlet", phi_sides);
SetNumber("Boundary Conditions/bottom_edge/dirichlet", phi_bottom);
SetNumber("Boundary Conditions/right_edge/dirichlet", phi_sides);
SetNumber("Boundary Conditions/top_edge/neumann", 0.);
\ No newline at end of file
BEM/l_shape.geo 0 → 100644
+ 55
− 0
View file @ d829d8ef
w_inlet = 0.5;
w_outlet = 1.0;
L_inlet = 1.0;
L_outlet = 1.5;
phi_inlet = 200.;
phi_outlet = 100.;
n_trans_inlet = 10;
n_long_inlet = 10;
n_trans_outlet = 10;
n_long_outlet = 15;
Point(1) = {0, 0, 0, 1.0};
Point(2) = {w_inlet, 0, 0, 1.0};
Point(3) = {w_inlet + L_outlet, 0, 0, 1.0};
Point(4) = {w_inlet + L_outlet, w_outlet, 0, 1.0};
Point(5) = {w_inlet, w_outlet, 0, 1.0};
Point(6) = {w_inlet, w_outlet + L_inlet, 0, 1.0};
Point(7) = {0, w_outlet + L_inlet, 0, 1.0};
Point(8) = {0, w_outlet, 0, 1.0};
Line(1) = {1, 2};
Line(2) = {2, 3};
Line(3) = {3, 4};
Line(4) = {4, 5};
Line(5) = {5, 6};
Line(6) = {6, 7};
Line(7) = {7, 8};
Line(8) = {8, 1};
Curve Loop(1) = {8, 1, 2, 3, 4, 5, 6, 7};
Plane Surface(1) = {1};
Transfinite Curve {1} = n_trans_inlet + 1 Using Progression 1;
Transfinite Curve {2} = n_long_outlet + 1 Using Progression 1;
Transfinite Curve {3} = n_trans_outlet + 1 Using Progression 1;
Transfinite Curve {4} = n_long_outlet + 1 Using Progression 1;
Transfinite Curve {5} = n_long_inlet + 1 Using Progression 1;
Transfinite Curve {6} = n_trans_inlet + 1 Using Progression 1;
Transfinite Curve {7} = n_long_inlet + 1 Using Progression 1;
Transfinite Curve {8} = n_trans_outlet + 1 Using Progression 1;
Recombine Surface {1};
Mesh.ElementOrder = 2;
Physical Curve("sides", 5) = {1, 2, 4, 5, 7, 8};
Physical Curve("inlet", 6) = {3};
Physical Curve("outlet", 7) = {6};
SetNumber("Boundary Conditions/sides/neumann", 0.);
SetNumber("Boundary Conditions/inlet/dirichlet", phi_inlet);
SetNumber("Boundary Conditions/outlet/dirichlet", phi_outlet);
\ No newline at end of file
BEM/main.cpp
+ 54
− 10
View file @ d829d8ef
//
// main.cpp
// BEM
//
// Created by Maltez Cavalheiro Kévin on 03/04/2022.
//
#include <gmsh.h>
#include <iostream>
#include <map>
#include <sstream>
#include <Eigen/Dense>
#include <cmath>
#include <iomanip>
# include <stdio.h>
#include <chrono>
#include "functions.hpp"
#include <omp.h>
using namespace std;
using namespace gmsh;
using namespace Eigen;
int main(int argc, char **argv)
{
if (argc < 2)
{
cout << "Usage: " << argv[0] << " <geo_file>\n";
return 0;
}
initialize();
open(argv[1]);
model::mesh::generate(2);
// #pragma omp parallel
// #pragma omp critical
// cout << "Greetings from thread "<< omp_get_thread_num() << endl;
vector<int> elementTypes;
vector<std::vector<size_t>> elementTags;
vector<std::vector<size_t>> nodeTags;
model::mesh::getElements(elementTypes, elementTags, nodeTags, 1);
int nb1dElements = 0;
for(int i = 0; i < elementTypes.size(); i++)
{
nb1dElements += elementTags[i].size();
}
vector<size_t> finalElementTags(nb1dElements);
vector<double> electrostaticPressure(nb1dElements);
solverBEMFun(argc, argv, finalElementTags, electrostaticPressure);
// for(int i = 0; i < finalElementTags.size(); i++)
// {
// cout << "Element tag '" << finalElementTags[i] << "'\t-> Electrostatic pressure = " << electrostaticPressure[i] << endl;
// }
finalize();
int main(int argc, const char * argv[]) {
// insert code here...
std::cout << "Hello, World!\n";
return 0;
}
BEM/my_geo.geo 0 → 100644
+ 51
− 0
View file @ d829d8ef
Lx = 1;
Ly = 1;
nx = 100;
ny = 100;
Point(1) = {0, 0, 0, 1.0};
Point(2) = {Lx, 0, 0, 1.0};
Point(3) = {Lx, Ly, 0, 1.0};
Point(4) = {0, Ly, 0, 1.0};
Line(1) = {1, 2};
Line(2) = {2, 3};
Line(3) = {3, 4};
Line(4) = {4, 1};
Curve Loop(1) = {4, 1, 2, 3};
Plane Surface(1) = {1};
Transfinite Curve {3, 1} = nx+1 Using Progression 1;
Transfinite Curve {4, 2} = ny+1 Using Progression 1;
Transfinite Surface {1};
Recombine Surface {1}; // quads instead of triangles
//à gauche c'est le tag qu'on donne, il doit être différent pour chaque curve ou surface et à droite c'est l'élément définit correspondant 4 courbes et 1 surface
Physical Curve("bottom_edge", 1) = {1};
Physical Curve("right_edge", 2) = {2};
Physical Curve("top_edge", 3) = {3};
Physical Curve("left_edge", 4) = {4};
//Physical Curve("boundary", 5) = {1,2,3,4};
Physical Surface("domain", 6) = {1};
// additional parameters given to the solver
// BC domain
//SetNumber("Materials/domain/Young", 210e3);
//SetNumber("Materials/domain/Poisson", 0.3);
//SetNumber("Materials/domain/rho",7800); //volumic mass of acier
//SetNumber("Volumic Forces/domain/bx",0.);
//SetNumber("Volumic Forces/domain/by",0.); //set to -9.81 for gravity
//BC top edge
SetNumber("Boundary Conditions/top_edge/neumann", 0.);
//BC right edge
SetNumber("Boundary Conditions/right_edge/dirichlet", 200);
//BC bottom edge
SetNumber("Boundary Conditions/bottom_edge/neumann", 0.);
// BC left edge
SetNumber("Boundary Conditions/left_edge/dirichlet", 100);
//+
Show "*";
BEM/rectangle.geo 0 → 100644
+ 35
− 0
View file @ d829d8ef
Lx = 1.5;
Ly = 1.0;
phi_left = 200.;
phi_right = 100.;
nx = 20;
ny = 5;
Point(1) = {0, 0, 0, 1.0};
Point(2) = {Lx, 0, 0, 1.0};
Point(3) = {Lx, Ly, 0, 1.0};
Point(4) = {0, Ly, 0, 1.0};
Line(1) = {1, 2};
Line(2) = {2, 3};
Line(3) = {3, 4};
Line(4) = {4, 1};
Curve Loop(1) = {4, 1, 2, 3};
Plane Surface(1) = {1};
Transfinite Curve {3, 1} = nx+1 Using Progression 1;
Transfinite Curve {4, 2} = ny+1 Using Progression 1;
Transfinite Surface {1};
Recombine Surface {1}; // quads instead of triangles
Physical Curve("left_edge", 5) = {4};
Physical Curve("bottom_edge", 6) = {1};
Physical Curve("right_edge", 7) = {2};
Physical Curve("top_edge", 8) = {3};
// additional parameters given to the solver
SetNumber("Boundary Conditions/left_edge/dirichlet", phi_left);
SetNumber("Boundary Conditions/bottom_edge/neumann", 0.);
SetNumber("Boundary Conditions/right_edge/dirichlet", phi_right);
SetNumber("Boundary Conditions/top_edge/neumann", 0.);
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment