Commit 04ad9198 authored by thomaskummer's avatar thomaskummer
Browse files

intrepid lifev update

parent 726d3f47
......@@ -53,47 +53,8 @@ using namespace LifeV;
// Functions
//============================================//
Real
externalPower ( const VectorEpetra& dispCurrent,
const VectorEpetra& dispPrevious,
const boost::shared_ptr<ETFESpace< RegionMesh<LinearTetra>, MapEpetra, 3, 3 > > dispETFESpace,
Real pressure,
Real dt,
const unsigned int bdFlag)
{
VectorEpetra traction ( dispCurrent.map() );
VectorEpetra velocity ( (dispCurrent - dispPrevious) / dt );
MatrixSmall<3,3> Id;
Id(0,0) = 1.; Id(0,1) = 0., Id(0,2) = 0.;
Id(1,0) = 0.; Id(1,1) = 1., Id(1,2) = 0.;
Id(2,0) = 0.; Id(2,1) = 0., Id(2,2) = 1.;
{
using namespace ExpressionAssembly;
auto I = value(Id);
auto Grad_u = grad( dispETFESpace, dispCurrent, 0);
auto F = Grad_u + I;
auto FmT = minusT(F);
auto J = det(F);
auto p = value(pressure);
QuadratureBoundary myBDQR (buildTetraBDQR (quadRuleTria7pt) );
integrate ( boundary ( dispETFESpace->mesh(), bdFlag),
myBDQR,
dispETFESpace,
p * J * dot( FmT * Nface, phi_i)
//p * J * dot( FmT * Nface, phi_i)
//value(-1.0) * J * dot (vE1, FmT * Nface) * phi_i) >> intergral
) >> traction;
traction.globalAssemble();
}
return traction.dot(velocity);
}
Real patchForce (const Real& t, const Real& Tmax, const Real& tmax, const Real& tduration)
{
......
......@@ -18,12 +18,13 @@
namespace LifeV
{
template<class meshType>
class HeartSolver {
public:
HeartSolver(EMSolver& emSolver, CirculationSolver& circulationSolver):
HeartSolver(EMSolver<meshType,EMMonodomainSolver<meshType> > & emSolver, CirculationSolver & circulationSolver):
M_emSolver (emSolver),
M_circulationSolver (circulationSolver)
{};
......@@ -37,7 +38,7 @@ public:
private:
protected:
EMSolver M_emSolver;
......@@ -50,6 +51,50 @@ private:
VectorSmall<2> M_volume;
Real
externalPower ( const VectorEpetra& dispCurrent,
const VectorEpetra& dispPrevious,
const boost::shared_ptr<ETFESpace< RegionMesh<LinearTetra>, MapEpetra, 3, 3 > > dispETFESpace,
Real pressure,
Real dt,
const unsigned int bdFlag) const
{
VectorEpetra traction ( dispCurrent.map() );
VectorEpetra velocity ( (dispCurrent - dispPrevious) / dt );
MatrixSmall<3,3> Id;
Id(0,0) = 1.; Id(0,1) = 0., Id(0,2) = 0.;
Id(1,0) = 0.; Id(1,1) = 1., Id(1,2) = 0.;
Id(2,0) = 0.; Id(2,1) = 0., Id(2,2) = 1.;
{
using namespace ExpressionAssembly;
auto I = value(Id);
auto Grad_u = grad( dispETFESpace, dispCurrent, 0);
auto F = Grad_u + I;
auto FmT = minusT(F);
auto J = det(F);
auto p = value(pressure);
QuadratureBoundary myBDQR (buildTetraBDQR (quadRuleTria7pt) );
integrate ( boundary ( dispETFESpace->mesh(), bdFlag),
myBDQR,
dispETFESpace,
p * J * dot( FmT * Nface, phi_i)
//p * J * dot( FmT * Nface, phi_i)
//value(-1.0) * J * dot (vE1, FmT * Nface) * phi_i) >> intergral
) >> traction;
traction.globalAssemble();
}
return traction.dot(velocity);
}
};
}
......
Supports Markdown
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