#include <SyFi.h>
#include <fstream>
Go to the source code of this file.
Functions | |
void | compute_poisson_element_matrix (FE &fe, Dof &dof, std::map< std::pair< unsigned int, unsigned int >, ex > &A) |
void | compute_nlconvdiff_element_matrix (FE &fe, Dof &dof, std::map< std::pair< unsigned int, unsigned int >, ex > &A) |
int | main () |
void compute_nlconvdiff_element_matrix | ( | FE & | fe, | |
Dof & | dof, | |||
std::map< std::pair< unsigned int, unsigned int >, ex > & | A | |||
) |
Definition at line 36 of file nljacobian_ex.cpp.
References demos::poisson1::Aij, SyFi::FE::dof(), fem_sympy::Fi, SyFi::FE::get_polygon(), SyFi::Dof::glob_dof(), SyFi::grad(), SyFi::inner(), SyFi::Dof::insert_dof(), SyFi::Polygon::integrate(), SyFi::FE::N(), SyFi::FE::nbf(), and fem_sympy::uj.
Referenced by main().
00040 { 00041 std::pair<unsigned int,unsigned int> index; 00042 Polygon& domain = fe.get_polygon(); 00043 00044 // insert the local dofs into the global Dof object 00045 for (unsigned int i=0; i< fe.nbf() ; i++) { 00046 dof.insert_dof(1,i,fe.dof(i)); 00047 } 00048 00049 // create the local U field: U = sum_k u_k N_k 00050 ex UU = matrix(2,1,lst(0,0)); 00051 ex ujs = symbolic_matrix(1,fe.nbf(), "u"); 00052 for (unsigned int k=0; k< fe.nbf(); k++) { 00053 UU +=ujs.op(k)*fe.N(k); // U += u_k N_k 00054 } 00055 00056 //Get U represented as a matrix 00057 matrix U = ex_to<matrix>(UU.evalm()); 00058 00059 for (unsigned int i=0; i< fe.nbf() ; i++) { 00060 index.first = dof.glob_dof(fe.dof(i)); // fetch global dof associated with i 00061 00062 // First: the diffusion term in Fi 00063 ex gradU = grad(U); // compute the gradient of U 00064 ex Fi_diffusion = inner(gradU, grad(fe.N(i))); // inner product of grad(U) and grad(Ni) 00065 00066 // Second: the convection term in Fi 00067 ex Ut = U.transpose(); // get the transposed of U 00068 ex UgradU = (Ut*gradU).evalm(); // compute U*grad(U) 00069 ex Fi_convection = inner(UgradU, fe.N(i), true); // compute U*grad(U)*Ni 00070 00071 // add together terms for convection and diffusion 00072 ex Fi = Fi_convection + Fi_diffusion; 00073 00074 00075 // Loop over all uj and differentiate Fi with respect 00076 // to uj to get the Jacobian Jij 00077 for (unsigned int j=0; j< fe.nbf() ; j++) { 00078 index.second = dof.glob_dof(fe.dof(j)); // fetch global dof associated with j 00079 symbol uj = ex_to<symbol>(ujs.op(j)); // cast uj to a symbol 00080 ex Jij = Fi.diff(uj,1); // differentiate Fi with respect to uj 00081 ex Aij = domain.integrate(Jij); // intergrate the Jacobian Jij 00082 A[index] += Aij; // update the global matrix 00083 } 00084 } 00085 }
void compute_poisson_element_matrix | ( | FE & | fe, | |
Dof & | dof, | |||
std::map< std::pair< unsigned int, unsigned int >, ex > & | A | |||
) |
Definition at line 8 of file nljacobian_ex.cpp.
References demos::poisson1::Aij, fem_sympy::Fi, SyFi::FE::get_polygon(), SyFi::grad(), SyFi::inner(), SyFi::Polygon::integrate(), SyFi::FE::N(), SyFi::FE::nbf(), fem_sympy::u, and fem_sympy::uj.
Referenced by main().
00012 { 00013 std::pair<unsigned int,unsigned int> index; 00014 Polygon& domain = fe.get_polygon(); 00015 00016 ex ujs = symbolic_matrix(1,fe.nbf(), "u"); 00017 ex u; 00018 for (unsigned int k=0; k< fe.nbf(); k++) { 00019 u += ujs.op(k)*fe.N(k); 00020 } 00021 00022 00023 for (unsigned int i=0; i< fe.nbf() ; i++) { 00024 index.first = i; 00025 ex Fi = inner(grad(u), grad(fe.N(i))); 00026 for (unsigned int j=0; j< fe.nbf() ; j++) { 00027 index.second = j; 00028 symbol uj = ex_to<symbol>(ujs.op(j)); 00029 ex nabla = Fi.diff(uj,1); 00030 ex Aij = domain.integrate(nabla); 00031 A[index] += Aij; 00032 } 00033 } 00034 }
int main | ( | ) |
Definition at line 90 of file nljacobian_ex.cpp.
References SyFi::compare_archives(), SyFi::VectorLagrange::compute_basis_functions(), SyFi::Lagrange::compute_basis_functions(), compute_nlconvdiff_element_matrix(), compute_poisson_element_matrix(), demos::crouzeixraviart::fe, SyFi::initSyFi(), SyFi::istr(), print(), SyFi::StandardFE::set_order(), SyFi::StandardFE::set_polygon(), SyFi::VectorLagrange::set_size(), and test::usage.
00090 { 00091 00092 initSyFi(2); 00093 00094 Triangle T(lst(0,0), lst(1,0), lst(0,1), "t"); 00095 int order = 2; 00096 00097 // First we compute a standard Poisson problem, i.e., 00098 // The differentiation of F(u_i) with respect to u_j 00099 // should give the standard Poisson problem. 00100 Lagrange fe; 00101 fe.set_order(order); 00102 fe.set_polygon(T); 00103 fe.compute_basis_functions(); 00104 00105 Dof dof1; 00106 std::map<std::pair<unsigned int,unsigned int>, ex> A1; 00107 compute_poisson_element_matrix(fe,dof1,A1); 00108 print(A1); 00109 00110 // Second we compute a nonlinear convection 00111 // diffusion problem. 00112 VectorLagrange vfe; 00113 vfe.set_order(order); 00114 vfe.set_size(2); 00115 vfe.set_polygon(T); 00116 vfe.compute_basis_functions(); 00117 usage(vfe); 00118 00119 Dof dof2; 00120 std::map<std::pair<unsigned int,unsigned int>, ex> A2; 00121 compute_nlconvdiff_element_matrix(vfe,dof2,A2); 00122 cout <<"standard output"<<endl; 00123 print(A2); 00124 cout <<"LaTeX output"<<endl; 00125 cout <<latex; 00126 print(A2); 00127 cout <<"Python output"<<endl; 00128 cout <<python; 00129 print(A2); 00130 cout <<"C output"<<endl; 00131 cout <<csrc; 00132 print(A2); 00133 00134 00135 00136 // regression test 00137 00138 archive ar; 00139 map<std::pair<unsigned int,unsigned int>,ex>::iterator iter; 00140 for (iter = A1.begin(); iter != A1.end() ; iter++) { 00141 ar.archive_ex((*iter).second, istr("A1_", 00142 (*iter).first.first, 00143 (*iter).first.second).c_str()); 00144 } 00145 for (iter = A2.begin(); iter != A2.end() ; iter++) { 00146 ar.archive_ex((*iter).second, istr("A2_", 00147 (*iter).first.first, 00148 (*iter).first.second).c_str()); 00149 } 00150 00151 ofstream vfile("nljacobian_ex.gar.v"); 00152 vfile << ar; vfile.close(); 00153 if(!compare_archives("nljacobian_ex.gar.v", "nljacobian_ex.gar.r")) { 00154 cerr << "Failure!" << endl; 00155 return -1; 00156 } 00157 00158 return 0; 00159 }