Trilinos based (stochastic) FEM solvers
main.cpp
Go to the documentation of this file.
1 /*
2 Brian Staber (brian.staber@gmail.com)
3 */
4 
5 #include "Epetra_ConfigDefs.h"
6 #ifdef HAVE_MPI
7 #include "mpi.h"
8 #include "Epetra_MpiComm.h"
9 #else
10 #include "Epetra_SerialComm.h"
11 #endif
12 
14 #include "Teuchos_RCP.hpp"
15 #include "Teuchos_ParameterList.hpp"
16 #include "Teuchos_StandardCatchMacros.hpp"
17 #include "Teuchos_CommandLineProcessor.hpp"
18 #include "Teuchos_XMLParameterListCoreHelpers.hpp"
19 
20 int main(int argc, char *argv[]){
21 
22  std::string xmlInFileName = "";
23 
24  Teuchos::CommandLineProcessor clp(false);
25  clp.setOption("xml-in-file",&xmlInFileName,"The XML file to read into a parameter list");
26  clp.setDocString("TO DO.");
27 
28  Teuchos::CommandLineProcessor::EParseCommandLineReturn
29  parse_return = clp.parse(argc,argv);
30  if( parse_return != Teuchos::CommandLineProcessor::PARSE_SUCCESSFUL ) {
31  std::cout << "\nEnd Result: TEST FAILED" << std::endl;
32  return parse_return;
33  }
34 
35 #ifdef HAVE_MPI
36 MPI_Init(&argc, &argv);
37  Epetra_MpiComm Comm(MPI_COMM_WORLD);
38 #else
39  Epetra_SerialComm Comm;
40 #endif
41 
42  Teuchos::RCP<Teuchos::ParameterList> paramList = Teuchos::rcp(new Teuchos::ParameterList);
43  if(xmlInFileName.length()) {
44  Teuchos::updateParametersFromXmlFile(xmlInFileName, inoutArg(*paramList));
45  }
46 
47  if (Comm.MyPID()==0){
48  paramList->print(std::cout,2,true,true);
49  }
50 
51  Epetra_SerialDenseVector x(7);
52  x(0) = Teuchos::getParameter<double>(paramList->sublist("TIMooney"),"mu1");
53  x(1) = Teuchos::getParameter<double>(paramList->sublist("TIMooney"),"mu2");
54  x(2) = Teuchos::getParameter<double>(paramList->sublist("TIMooney"),"mu3");
55  x(3) = Teuchos::getParameter<double>(paramList->sublist("TIMooney"),"mu4");
56  x(4) = Teuchos::getParameter<double>(paramList->sublist("TIMooney"),"mu5");
57  x(5) = Teuchos::getParameter<double>(paramList->sublist("TIMooney"),"beta4");
58  x(6) = Teuchos::getParameter<double>(paramList->sublist("TIMooney"),"beta5");
59  for (unsigned int i=0; i<5; i++){
60  x(i) = 1.0e3*x(i);
61  }
62  double angle = 15.0*2.0*M_PI/360.0;
63 
64  Teuchos::RCP<compressibleHyperelasticity_linearPatchTest> obj =
65  Teuchos::rcp(new compressibleHyperelasticity_linearPatchTest(Comm,*paramList));
66  obj->set_parameters(x,angle);
67 
68  Teuchos::RCP<newtonRaphson> newton =
69  Teuchos::rcp(new newtonRaphson(*obj,*paramList));
70 
71  newton->Initialization();
72  int flagNewton = newton->Solve_with_Aztec(true);
73 
74  double errorL2 = 1.0;
75  if (!flagNewton){
76  errorL2 = obj->errorL2(*newton->x);
77  std::string filenameU = "/Users/brian/Documents/GitHub/Trilinos_results/nrl/UlinearPatchTest.mtx";
78  newton->print_newton_solution(filenameU);
79  if(Comm.MyPID()==0){
80  std::cout << "-----------------------------\n";
81  std::cout << "||u||_L^2 = " << errorL2 << "\n";
82  std::cout << "-----------------------------\n";
83  }
84  }
85 
86 #ifdef HAVE_MPI
87  MPI_Finalize();
88 #endif
89 return 0;
90 
91 }
int main(int argc, char *argv[])
Definition: main.cpp:23
for i
Definition: costFunction.m:38