00001 #include "biotfem.h"
00002 #include "assemblyorthomeshsystem.h"
00003 #include "hdf5orthowriter.h"
00004 #include "gnuplotanim.h"
00005
00006
00007 BiotFEM::BiotFEM(OrthoMesh &mesh, double dt, VecDouble &cK,double &fYoung, double &fPoisson,Function3D &fD, Function3D &fN)
00008 :m_mesh(mesh),m_cK(cK),m_eqvar(mesh,cK,fYoung,fPoisson,fD,fN),m_fN(fN)
00009 {
00010 m_eqvar.printDealFEInfo();
00011 exit(0);
00012 AssemblyOrthoMeshSystem<IsotropicElasticEqVar> assembler;
00013 m_vRHS.reinit(m_eqvar.n_dofs());
00014 m_sol.reinit(m_eqvar.n_dofs());
00015 m_sparsity.reinit(m_eqvar.n_dofs(),m_eqvar.n_dofs(),m_eqvar.n_max_dof_coupling());
00016 assembler.buildSparsityPattern(m_sparsity,m_mesh,m_eqvar);
00017 m_sparsity.compress();
00018 m_M.reinit(m_sparsity);
00019 assembler.setBoundaryValue(m_mesh,m_eqvar,fD,dirichlet_bc_map);
00020
00021 }
00022
00023
00024
00025 void BiotFEM::iterate(TransportBase &trans)
00026 {
00027 AssemblyOrthoMeshSystem<IsotropicElasticEqVar> assembler;
00028 assembler.buildGlobalMatrix(m_M,m_mesh, m_eqvar);
00029 assembler.buildGlobalRhs(m_vRHS,m_mesh,m_eqvar);
00030 assembler.applyBoundaryValue(dirichlet_bc_map,m_M);
00031 assembler.applyBoundaryValue(dirichlet_bc_map,m_vRHS);
00032
00033 m_sol=m_vRHS;
00034 m_sol.print();
00035 solveA.factorize(m_M);
00036 solveA.solve(m_sol);
00037
00038 }
00039 void BiotFEM::printOutput()
00040 {
00041 return;
00042 if (m_debugLevel > 0)
00043 {
00044 HDF5OrthoWriter &hdf5 = HDF5OrthoWriter::getHDF5OrthoWriter();
00045 hdf5.writeScalarField(m_sol,"P");
00046
00047 GnuPlotAnim &plot = GnuPlotAnim::getGnuPlotAnim();
00048
00049
00050
00051 plot.set("yrange","[0:1]");
00052 plot.plotVerticeValuesAtFixedZ(m_mesh,m_sol,"Sc","with lines",5);
00053 plot.NextScene();
00054 }
00055 }
00056
00057
00058 void BiotFEM::getVelocitiesAtFaces(Matrix& M)
00059 {
00060 throw new Exception("BiotFEM::getVelocitiesAtFaces Not implemented for this class\n");
00061 }