CArl
Code Arlequin / C++ implementation
libmesh_assemble_lin_homogeneous.cpp File Reference

Go to the source code of this file.

Functions

int main (int argc, char **argv)
 

Function Documentation

int main ( int  argc,
char **  argv 
)

Definition at line 3 of file libmesh_assemble_lin_homogeneous.cpp.

3  {
4 
5  // [USER] Fixed boudary
7 
8  // --- Initialize libMesh
9  libMesh::LibMeshInit init(argc, argv);
10 
11  // Do performance log?
12  const bool MASTER_bPerfLog_carl_libmesh = true;
13  libMesh::PerfLog perf_log("Main program", MASTER_bPerfLog_carl_libmesh);
14 
15  // libMesh's C++ / MPI communicator wrapper
16  libMesh::Parallel::Communicator& WorldComm = init.comm();
17 
18  // Number of processors and processor rank.
19  int rank = WorldComm.rank();
20  int nodes = WorldComm.size();
21 
22  // --- Set up inputs
23 
24  // Command line parser
25  GetPot command_line(argc, argv);
26 
27  // File parser
28  GetPot field_parser;
29 
30  // If there is an input file, parse it to get the parameters. Else, parse the command line
31  std::string input_filename;
32  if (command_line.search(2, "--inputfile", "-i")) {
33  input_filename = command_line.next(input_filename);
34  field_parser.parse_input_file(input_filename, "#", "\n", " \t\n");
35  } else {
36  field_parser = command_line;
37  }
38 
39  libmesh_assemble_input_params input_params;
40  get_input_params(field_parser, input_params);
41 
42  // Check libMesh installation dimension
43  const unsigned int dim = 3;
44 
45  libmesh_example_requires(dim == LIBMESH_DIM, "3D support");
46 
47  // --- Declare the three meshes to be intersected
48 
49  // - Parallelized meshes: A, B, mediator and weight
50  perf_log.push("Meshes - Parallel","Read files:");
51  libMesh::Mesh system_mesh(WorldComm, dim);
52  system_mesh.read(input_params.mesh_file);
53  system_mesh.prepare_for_use();
54 
55  libMesh::Mesh mesh_weight(WorldComm, dim);
56  mesh_weight.allow_renumbering(false);
57  mesh_weight.read(input_params.mesh_weight_file);
58  mesh_weight.prepare_for_use();
59 
60  perf_log.pop("Meshes - Parallel","Read files:");
61 
62  // --- Generate the equation systems
63  perf_log.push("System setup:");
64 
65  libMesh::EquationSystems& equation_systems = libMesh::EquationSystems(system_mesh);
66 
67  // Add linear elasticity system
68  libMesh::LinearImplicitSystem& elasticity_system
69  = add_elasticity(equation_systems);
70 
71  // Set clamped border
72  set_clamped_border(elasticity_system, fixed_bound_id);
73 
74 // // [MACRO] Build stress system
75 // libMesh::ExplicitSystem& stress_system_BIG
76 // = add_stress(equation_systems_BIG);
77 
78 // equation_systems_BIG.init();
79 
80 // perf_log.pop("Macro system","System initialization:");
81 
82 // // - Build the micro system
83 
84 // perf_log.push("Micro system","System initialization:");
85 // libMesh::EquationSystems& equation_systems_micro =
86 // CoupledTest.add_micro_EquationSystem<libMesh::PetscMatrix<libMesh::Number> >("MicroSys", mesh_micro);
87 
88 // // [MICRO] Set up the physical properties
89 // libMesh::LinearImplicitSystem& elasticity_system_micro
90 // = add_elasticity(equation_systems_micro);
91 
92 // // [MICRO] Build stress system
93 // libMesh::ExplicitSystem& stress_system_micro
94 // = add_stress(equation_systems_micro);
95 
96 // equation_systems_micro.init();
97 // perf_log.pop("Micro system","System initialization:");
98 
99 // // - Build the RESTRICTED BIG system
100 // perf_log.push("RESTRICTED macro system","System initialization:");
101 // libMesh::EquationSystems& equation_systems_R_BIG =
102 // CoupledTest.set_Restricted_BIG_EquationSystem("BigSys", mesh_R_BIG);
103 
104 // // [R. MACRO] Set up the physical properties
105 // libMesh::ExplicitSystem& elasticity_system_R_BIG
106 // = add_explicit_elasticity(equation_systems_R_BIG);
107 
108 // carl::reduced_system_init(elasticity_system_R_BIG);
109 
110 // perf_log.pop("RESTRICTED macro system","System initialization:");
111 
112 // // - Build the RESTRICTED micro system
113 
114 // perf_log.push("RESTRICTED micro system","System initialization:");
115 // libMesh::EquationSystems& equation_systems_R_micro =
116 // CoupledTest.add_Restricted_micro_EquationSystem("MicroSys", mesh_R_micro);
117 
118 // // [R. MICRO] Set up the physical properties
119 // libMesh::ExplicitSystem& elasticity_system_R_micro
120 // = add_explicit_elasticity(equation_systems_R_micro);
121 
122 // carl::reduced_system_init(elasticity_system_R_micro);
123 // perf_log.pop("RESTRICTED micro system","System initialization:");
124 
125 // // - Build the mediator system
126 
127 // perf_log.push("Mediator system","System initialization:");
128 // libMesh::EquationSystems& equation_systems_mediator =
129 // CoupledTest.add_mediator_EquationSystem("MediatorSys", mesh_mediator);
130 
131 // libMesh::LinearImplicitSystem& elasticity_system_mediator
132 // = add_elasticity(equation_systems_mediator);
133 
134 // equation_systems_mediator.init();
135 
136 // WorldComm.barrier();
137 
138 // perf_log.pop("Mediator system","System initialization:");
139 
140 // // - Build the dummy inter system
141 
142 // perf_log.push("Intersection system","System initialization:");
143 // libMesh::ExplicitSystem& elasticity_system_inter
144 // = add_explicit_elasticity(equation_systems_inter);
145 
146 // carl::reduced_system_init(elasticity_system_inter);
147 
148 // WorldComm.barrier();
149 
150 // perf_log.pop("Intersection system","System initialization:");
151 
152 // perf_log.push("Physical properties","System initialization:");
153 
154 // // --- Set the physical parameters
155 // double BIG_E = 0;
156 // double BIG_Mu = 0;
157 
158 // double coupling_const = -1;
159 
160 // // Anisotropic properties for the micro system
161 // carl::anisotropic_elasticity_tensor_cubic_sym anisotropy_data(equation_systems_micro,input_params.physical_params_file,BIG_E,BIG_Mu);
162 
163 // // Homogeneous properties for the macro system
164 // set_constant_physical_properties(equation_systems_BIG,BIG_E,BIG_Mu);
165 
166 // perf_log.pop("Physical properties","System initialization:");
167 
168 // // --- Set the coupling matrices
169 // perf_log.push("Set coupling matrices");
170 // coupling_const = BIG_E;
171 
172 // // Set parameters
173 // CoupledTest.set_coupling_parameters("MicroSys",coupling_const,input_params.mean_distance);
174 
175 // // Use H1 norm
176 // CoupledTest.use_H1_coupling("MicroSys");
177 
178 // // Assemble!
179 // CoupledTest.assemble_coupling_elasticity_3D_parallel("BigSys","MicroSys",
180 // "InterSys","MediatorSys",
181 // mesh_R_BIG, mesh_R_micro,
182 // local_intersection_pairs_map,
183 // local_intersection_restricted_pairs_map,
184 // local_intersection_meshI_to_inter_map,
185 // inter_mediator_BIG,
186 // inter_mediator_micro);
187 
188 // // Print some info
189 // std::cout << std::endl;
190 // std::cout << "| ---> Constants " << std::endl;
191 // std::cout << "| Macro :" << std::endl;
192 // std::cout << "| E : " << BIG_E << std::endl;
193 // std::cout << "| Mu (lamba_2) : " << BIG_Mu << std::endl;
194 // std::cout << "| lambda_1 : " << eval_lambda_1(BIG_E,BIG_Mu) << std::endl;
195 // std::cout << "| Coupling :" << std::endl;
196 // std::cout << "| kappa : " << coupling_const << std::endl;
197 // std::cout << "| e : " << input_params.mean_distance << std::endl;
198 
199 // switch(input_params.solver_type)
200 // {
201 // case carl::LATIN_MODIFIED_STIFFNESS:
202 // case carl::LATIN_ORIGINAL_STIFFNESS:
203 // {
204 // std::cout << "| LATIN :" << std::endl;
205 // std::cout << "| k_dA, k_dB : " << input_params.k_dA << " " << input_params.k_dB << std::endl;
206 // std::cout << "| k_cA, k_cB : " << input_params.k_cA << " " << input_params.k_cB << std::endl;
207 // break;
208 // }
209 // case carl::CG:
210 // {
211 // break;
212 // }
213 // }
214 // perf_log.pop("Set coupling matrices");
215 
216 // std::cout << "| restart file : " << input_params.coupled_restart_file_base << "*" << std::endl;
217 
218 // std::cout << std::endl << "| --> Testing the solver " << std::endl << std::endl;
219 // perf_log.push("Set up","Coupled Solver:");
220 
221 // // --- Set restart files
222 // if(input_params.b_PrintRestartFiles || input_params.b_UseRestartFiles)
223 // {
224 // CoupledTest.set_restart( input_params.b_UseRestartFiles,
225 // input_params.b_PrintRestartFiles,
226 // input_params.coupled_restart_file_base);
227 // }
228 // std::cout << std::endl << "| --> Setting the solver " << std::endl << std::endl;
229 
230 // // --- Set up the coupled solver!
231 // switch(input_params.solver_type)
232 // {
233 // case carl::LATIN_MODIFIED_STIFFNESS:
234 // case carl::LATIN_ORIGINAL_STIFFNESS:
235 // {
236 // CoupledTest.set_LATIN_solver( "MicroSys","Elasticity",
237 // assemble_elasticity_with_weight,
238 // assemble_elasticity_anisotropic_with_weight,
239 // anisotropy_data,
240 // input_params.k_dA, input_params.k_dB, input_params.k_cA, input_params.k_cB,
241 // input_params.LATIN_eps, input_params.LATIN_conv_max, input_params.LATIN_relax);
242 // break;
243 // }
244 // case carl::CG:
245 // {
246 // CoupledTest.use_null_space_micro("MicroSys",true);
247 // CoupledTest.set_cg_preconditioner_type(input_params.CG_precond_type);
248 // CoupledTest.set_CG_solver( "MicroSys","Elasticity",
249 // assemble_elasticity_with_weight,
250 // assemble_elasticity_anisotropic_with_weight,
251 // anisotropy_data,
252 // input_params.CG_coupled_conv_abs,input_params.CG_coupled_conv_rel,
253 // input_params.CG_coupled_conv_max,input_params.CG_coupled_div,
254 // input_params.CG_coupled_conv_corr);
255 // break;
256 // }
257 // }
258 // perf_log.pop("Set up","Coupled Solver:");
259 
260 
261 // // --- Solve !
262 // perf_log.push("Solve","Coupled Solver:");
263 // CoupledTest.solve("MicroSys","Elasticity",input_params.coupled_convergence_output);
264 // perf_log.pop("Solve","Coupled Solver:");
265 
266 // // --- Calculate stress
267 // perf_log.push("Compute stress - micro","Output:");
268 // compute_stresses(equation_systems_micro);
269 // perf_log.pop("Compute stress - micro","Output:");
270 
271 // perf_log.push("Compute stress - macro","Output:");
272 // compute_stresses(equation_systems_BIG);
273 // perf_log.pop("Compute stress - macro","Output:");
274 
275 // // --- Export solution
276 // #ifdef LIBMESH_HAVE_EXODUS_API
277 // if(input_params.b_PrintOutput)
278 // {
279 // perf_log.push("Save output","Output:");
280 // libMesh::ExodusII_IO exo_io_micro(mesh_micro, /*single_precision=*/true);
281 
282 // std::set<std::string> system_names_micro;
283 // system_names_micro.insert("Elasticity");
284 // exo_io_micro.write_equation_systems(input_params.output_file_micro,equation_systems_micro,&system_names_micro);
285 
286 // exo_io_micro.write_element_data(equation_systems_micro);
287 
288 // libMesh::ExodusII_IO exo_io_BIG(system_mesh, /*single_precision=*/true);
289 
290 // std::set<std::string> system_names_BIG;
291 // system_names_BIG.insert("Elasticity");
292 // exo_io_BIG.write_equation_systems(input_params.output_file_BIG,equation_systems_BIG,&system_names_BIG);
293 
294 // exo_io_BIG.write_element_data(equation_systems_BIG);
295 // perf_log.pop("Save output","Output:");
296 // }
297 // #endif
298 
299 // if(input_params.b_ExportScalingData)
300 // {
301 // CoupledTest.print_perf_log(input_params.scaling_data_file);
302 // }
303 
304  return 0;
305 }
std::string mesh_weight_file
Path to the mesh containing the weight region indices.
void get_input_params(GetPot &field_parser, feti_iterate_params &input_params)
Parser function for the coupled solver test programs.
void set_clamped_border(libMesh::ImplicitSystem &elasticity_system, int boundary_id)
Set a clamped border.
libMesh::ImplicitSystem & add_elasticity(libMesh::EquationSystems &input_systems, libMesh::Order order=libMesh::FIRST, libMesh::FEFamily family=libMesh::LAGRANGE)
Add a linear elasticity libMesh::LinearImplicitSystem to the input libMesh::EquationSystems& input_sy...
boundary_id_cube
Small enumerate defining a cube's boundary ID's.
std::string mesh_file
Path to the system mesh.