SCoupling< M, O > Class Template Reference

#include <SCoupling.h>

Inheritance diagram for SCoupling< M, O >:

Classes
struct	SCoupledElement

Public Types
typedef O::ELEMENT	ELEMENT

Public Member Functions
	SCoupling ()=default
void	solveSurfaceCoupling (const unsigned max_adapt=0)
void	solveSurfaceCoupling (unsigned nStep, const unsigned max_adapt)
void	solveSurfaceCouplingForgiving (unsigned nStep, double timeMaxMin=-constants::inf, const unsigned max_adapt=0)
void	solveSurfaceCouplingFixedSolid ()
TriangleWall *	createTriangleWall (std::array< Vec3D, 3 > vertex)
void	updateTriangleWall (TriangleWall *&wall, std::array< Vec3D, 3 > vertex)
void	computeOneTimeStepForSCoupling (const unsigned &nStepsMercury, const unsigned max_adapt=0)
void	createDPMWallsFromFiniteElems ()
void	updateDPMWallsFromFiniteElems ()
void	updateTractionOnFiniteElems ()
bool	computeSCouplingForcesFromTriangles (ELEMENT *const elem_pt, const unsigned &nTriangles, unsigned &wallID, Vector< Vector< double > > &nodalCouplingForces)
void	getElementBoundingBox (ELEMENT *&elem_pt, Vec3D &min, Vec3D &max)
void	getSCoupledElements ()
void	coupleBoundary (unsigned b)
void	coupleBoundaries (std::vector< unsigned > b)
void	disableLogSurfaceCoupling ()
void	setSolidFeelsParticles (bool val)
bool	getSolidFeelsParticles () const
Public Member Functions inherited from BaseCoupling< M, O >
	BaseCoupling ()=default
void	setName (std::string name)
std::string	getName () const
void	removeOldFiles () const
void	writeEneTimeStep (std::ostream &os) const override
void	writeEneHeader (std::ostream &os) const override
void	solveOomph (int max_adapt=0)
void	checkResidual ()
void	solveMercury (unsigned long nt)
void	setCGWidth (const double &width)
double	getCGWidth ()
bool	useCGMapping ()
CGFunctions::LucyXYZ	getCGFunction ()

Private Attributes
Vector< SCoupledElement >	sCoupledElements_
	List of surface-coupled elements. More...
std::vector< TriangleWall * >	triangleWalls_
	List of triangle walls used for the surface coupling. More...
std::vector< unsigned >	coupledBoundaries_
bool	logSurfaceCoupling = true
bool	solidFeelsParticles_ = true

Member Typedef Documentation

ELEMENT

template<class M , class O >

typedef O::ELEMENT SCoupling< M, O >::ELEMENT

Constructor & Destructor Documentation

SCoupling()

template<class M , class O >

SCoupling< M, O >::SCoupling ( )

default

Member Function Documentation

computeOneTimeStepForSCoupling()

template<class M , class O >

void SCoupling< M, O >::computeOneTimeStepForSCoupling ( const unsigned &  nStepsMercury, const unsigned  max_adapt = 0  )

inline

Solve surface-coupled problem for one oomph time step (n mercury time steps)

    {
        auto t0 = std::chrono::system_clock::now();
        updateDPMWallsFromFiniteElems();
        auto t1 = std::chrono::system_clock::now();
        BaseCoupling<M,O>::solveMercury(nStepsMercury);
        auto t2 = std::chrono::system_clock::now();
        if (solidFeelsParticles_) {
            updateTractionOnFiniteElems();
        }
        auto t3 = std::chrono::system_clock::now();
        BaseCoupling<M,O>::solveOomph(max_adapt);
        auto t4 = std::chrono::system_clock::now();
        if (logSurfaceCoupling) logger(INFO, "time % Elapsed time: FEM->DEM %, DEM %, DEM->FEM %, FEM %", M::getTime(),
               (t1-t0).count(), (t2-t1).count(), (t3-t2).count(), (t4-t3).count());
        if(max_adapt>0)
        {
            getSCoupledElements();
            createDPMWallsFromFiniteElems();
        }
    }

References INFO, logger, BaseCoupling< M, O >::solveMercury(), and BaseCoupling< M, O >::solveOomph().

computeSCouplingForcesFromTriangles()

template<class M , class O >

bool SCoupling< M, O >::computeSCouplingForcesFromTriangles ( ELEMENT *const  elem_pt, const unsigned &  nTriangles, unsigned &  wallID, Vector< Vector< double > > &  nodalCouplingForces  )

inline

Computes nodal scoupling forces from triangles.

Parameters

[in]	elem_pt	pointer to the element for which forces are computed
[in]	nTriangles	number of triangles (equals number of nodes of the traction element)
[in,out]	wallID
[out]	nodalCouplingForces[l]	= sum_i f * psi(l,cp) for each shape function l, interaction i, evaluated at contact point cp

Returns

    {
        // whether the element is coupled to particles
        bool elemIsCoupled = false;
        
        // get number of nodes in the bulk element
        const unsigned nnode = elem_pt->nnode();
        // get dimension of the problem
        const unsigned dim = elem_pt->dim();
        // initialize the coupling force vector
        nodalCouplingForces.resize(nnode, Vector<double>(dim, 0.0));
        
        // get number of TriangleWalls created from an oomph face element
        unsigned n = 0;
        // loop over TriangleWalls
        while (n++ < nTriangles)
        {
            // get pointer to the wall
            TriangleWall* w = triangleWalls_[wallID++];
            
            // skip the wall if not interactions with it
            if (w->getInteractions().size() == 0) continue;
            
            // Set up memory for the shape/test functions
            Shape psi(nnode);
            
            // loop over interactions with the wall
            for (auto inter : w->getInteractions())
            {
                if (!inter->getForce().isZero())
                {
                    // scale contact point and forces from DEM to FEM units
                    Vec3D xc = inter->getContactPoint();
                    Vec3D fc = inter->getForce();
                    
                    Vector<double> x(3, 0.0), f(3, 0.0);
                    x[0] = xc.getX();
                    x[1] = xc.getY();
                    x[2] = xc.getZ();
                    f[0] = fc.getX();
                    f[1] = fc.getY();
                    f[2] = fc.getZ();
                    
                    // get the local coordinate s if xc is located in the finite element
                    Vector<double> s(3, 0.0);
                    GeomObject* geom_obj_pt = 0;
                    elem_pt->locate_zeta(x, geom_obj_pt, s);
                    
                    // Get shape/test fcts
                    elem_pt->shape(s, psi);
                    // Loop over the test functions
                    for (unsigned l = 0; l < nnode; l++)
                    {
                        //Loop over the force components
                        for (unsigned i = 0; i < dim; i++)
                        {
                            // add contribution to the nodal coupling force
                            nodalCouplingForces[l][i] += f[i] * psi(l);
                        }
                    }
                    // set the flag to true
                    elemIsCoupled = true;
                }
            }
        }
        return elemIsCoupled;
    }

References f(), Vec3D::getX(), Vec3D::getY(), Vec3D::getZ(), i, oomph::GeomObject::locate_zeta(), n, s, w, plotDoE::x, and MultiOpt::xc.

coupleBoundaries()

template<class M , class O >

void SCoupling< M, O >::coupleBoundaries ( std::vector< unsigned >  b )

inline

                                                 {
        coupledBoundaries_ = std::move(b);
    }

References b.

coupleBoundary()

template<class M , class O >

void SCoupling< M, O >::coupleBoundary ( unsigned  b )

inline

                                    {
        coupledBoundaries_ = {b};
    }

References b.

Referenced by ElementAnalysis::ElementAnalysis(), CoupledBeam::setupOomph(), and CoupledProblem::setupOomph().

createDPMWallsFromFiniteElems()

template<class M , class O >

void SCoupling< M, O >::createDPMWallsFromFiniteElems ( )

inline

    {
        // loop over bulk elements at boundary
        for (auto sCoupledElement : sCoupledElements_)
        {
            // loop over nodes in element at boundary
            Vector<Vector<double*> > position_pt = sCoupledElement.surface_node_position_pt;
            // reordering vertices of oomph face element (default = {0,1,3,2})
            swap(position_pt[2], position_pt[3]);
            
            // get global coordinate at the center
            Vec3D center;
            Vector<double> x(3, 0.0);
            sCoupledElement.bulk_elem_pt->interpolated_x(sCoupledElement.center_local, x);
            center.setX(x[0]);
            center.setY(x[1]);
            center.setZ(x[2]);
            
            // create TriangleWalls from oomph face element
            const unsigned nTriangles = position_pt.size();
            unsigned n = 0;
            while (n < nTriangles)
            {
                // get vertices of TriangleWall (multiply vertex position with the length scale of the O)
                std::array<Vec3D, 3> vertex;
                // one vertex at the center
                vertex[0] = center;
                // two vertices from the O<element,TIMESTEPPER>
                vertex[1] = Vec3D(*position_pt[0][0],*position_pt[0][1],*position_pt[0][2]);
                vertex[2] = Vec3D(*position_pt[1][0],*position_pt[1][1],*position_pt[1][2]);
                
                // create triangle facet
                TriangleWall* w = createTriangleWall(vertex);
                triangleWalls_.push_back(w);
                
                // rotate forward by one element
                rotate(position_pt.begin(), position_pt.begin() + 1, position_pt.end());
                n++;
            }
        }
    }

References n, Vec3D::setX(), Vec3D::setY(), Vec3D::setZ(), swap(), w, and plotDoE::x.

createTriangleWall()

template<class M , class O >

TriangleWall* SCoupling< M, O >::createTriangleWall ( std::array< Vec3D, 3 >  vertex )

inline

Create a triangle wall from a set of vertices

Parameters

vertex

Todo:

Why we need to set GroupID. I don't remember anymore.

Returns

    {
        TriangleWall wall;
        auto species = M::speciesHandler.getObject(0);
        wall.setSpecies(species);
        wall.setGroupId(100);
        wall.setVertices(vertex[0], vertex[1], vertex[2]);
        auto w = M::wallHandler.copyAndAddObject(wall);
        return w;
    }

References BaseObject::setGroupId(), BaseWall::setSpecies(), TriangleWall::setVertices(), and w.

disableLogSurfaceCoupling()

template<class M , class O >

void SCoupling< M, O >::disableLogSurfaceCoupling ( )

inline

                                     {
        logSurfaceCoupling = false;
    }

getElementBoundingBox()

template<class M , class O >

void SCoupling< M, O >::getElementBoundingBox ( ELEMENT *&  elem_pt, Vec3D &  min, Vec3D &  max  )

inline

needed in computeSCouplingForcesFromCG to decide which particles affect a given element

    {
        // three arrays that contain the x, y and z coordinates of the bulk element
        Vector<double> listOfCoordX;
        Vector<double> listOfCoordY;
        Vector<double> listOfCoordZ;
        
        // get the x, y and z coordinates of the bulk element
        const unsigned nnode = elem_pt->nnode();
        for (unsigned n = 0; n < nnode; n++)
        {
            listOfCoordX.push_back(elem_pt->node_pt(n)->x(0));
            listOfCoordY.push_back(elem_pt->node_pt(n)->x(1));
            listOfCoordZ.push_back(elem_pt->node_pt(n)->x(2));
        }
        
        // get the bounding box of the bulk element
        min.X = *min_element(listOfCoordX.begin(), listOfCoordX.end());
        min.Y = *min_element(listOfCoordY.begin(), listOfCoordY.end());
        min.Z = *min_element(listOfCoordZ.begin(), listOfCoordZ.end());
        max.X = *max_element(listOfCoordX.begin(), listOfCoordX.end());
        max.Y = *max_element(listOfCoordY.begin(), listOfCoordY.end());
        max.Z = *max_element(listOfCoordZ.begin(), listOfCoordZ.end());
        
        // extend the bounding box if construct mapping with coarse graining
        logger.assert_always(M::particleHandler.getLargestParticle(), "No particles detected");
        //todo does it make sense with -2R here?
        min -= Vec3D(1.0, 1.0, 1.0) * ( BaseCoupling<M,O>::getCGWidth() - 2 * M::particleHandler.getLargestParticle()->getRadius() );
        max += Vec3D(1.0, 1.0, 1.0) * ( BaseCoupling<M,O>::getCGWidth() - 2 * M::particleHandler.getLargestParticle()->getRadius() );
    }

References BaseCoupling< M, O >::getCGWidth(), logger, max, min, and n.

getSCoupledElements()

template<class M , class O >

void SCoupling< M, O >::getSCoupledElements ( )

inline

Get bulk elements along boundaries (SCoupling). The results are stored in sCoupledElements_.

    {
        // first clear bulk elements (for refineable problem)
        sCoupledElements_.clear();
 
        // loop over all boundaries
        for (unsigned b : coupledBoundaries_)
        {
            // number of bulk elements adjacent to boundary b
            unsigned n_element = this->solid_mesh_pt()->nboundary_element(b);
            // we only need to couple the elements on the upper boundary
            logger(INFO,"Coupling boundary %. which has % elements", b, n_element);
 
            // loop over the bulk elements adjacent to boundary b
            for (unsigned e = 0; e < n_element; e++)
            {
                // initialize the struct SCoupledElement
                SCoupledElement sCoupledElement;
                
                // get pointer to the bulk element that is adjacent to boundary b
                sCoupledElement.bulk_elem_pt = dynamic_cast<ELEMENT*>(this->solid_mesh_pt()->boundary_element_pt(b, e));
                
                // get the index of the face of element e along boundary b
                sCoupledElement.face_index = this->solid_mesh_pt()->face_index_at_boundary(b, e);
                
                // create temporary traction element
                SolidTractionElement<ELEMENT> trac_elem(sCoupledElement.bulk_elem_pt, sCoupledElement.face_index);
                
                // number of nodes on traction element
                unsigned n_node = trac_elem.nnode();
                
                // allocate space to store the local coordinates of the traction element's vertices
                sCoupledElement.surface_node_position_pt.reserve(n_node);
                
                // store the local coordinates of the traction element's center
                sCoupledElement.center_local.resize(3, 0.0);
                
                // assign addresses of nodal coordinates to the struct SCoupledElement
                for (unsigned n = 0; n < n_node; n++)
                {
                    // store pointer to solid nodes
                    sCoupledElement.node_pt.push_back(dynamic_cast<CoupledSolidNode*>(trac_elem.node_pt(n)));
                    
                    Vector<double> s(2, 0.0);
                    trac_elem.local_coordinate_of_node(n, s);
                    
                    Vector<double> s_bulk(3, 0.0);
                    trac_elem.get_local_coordinate_in_bulk(s, s_bulk);
                    
                    // pass the addresses of x coordinates at boundary to x_pt
                    Vector<double*> x_pt;
                    for (unsigned i = 0; i < 3; i++)
                    {
                        x_pt.push_back(&trac_elem.node_pt(n)->x(i));
                        sCoupledElement.center_local[i] += s_bulk[i];
                    }
                    /*
                                        // debug that the addresses of x coordinates are passed to surface_node_position_pt
                                        cout << "Address of x(1) of node on traction element\t " << &trac_elem.node_pt(n)->x(0) << endl;
                                        int nodeIDOfBulkElement = sCoupledElement.bulk_elem_pt->get_node_number(trac_elem.node_pt(n));
                                        cout << "Address of x(1) of the same node in bulk element "
                                             << &sCoupledElement.bulk_elem_pt->node_pt(nodeIDOfBulkElement)->x(0) << endl;
                                        cout << "Address of x(1) saved in Oomph interface object\t " << x_pt[0] << endl;
                    */
                    // save the addresses of Lagrange coordinates in each element at boundary
                    sCoupledElement.surface_node_position_pt.push_back(x_pt);
                }
                
                // local coordinate of the center of the face element
                for (unsigned i = 0; i < 3; i++) sCoupledElement.center_local[i] /= n_node;
                
                // save bulk elements and surface nodal position pointers
                sCoupledElements_.push_back(sCoupledElement);
            }
        }
    }

References b, SCoupling< M, O >::SCoupledElement::bulk_elem_pt, SCoupling< M, O >::SCoupledElement::center_local, e(), SCoupling< M, O >::SCoupledElement::face_index, oomph::FaceElement::get_local_coordinate_in_bulk(), i, INFO, oomph::FiniteElement::local_coordinate_of_node(), logger, n, oomph::FiniteElement::nnode(), SCoupling< M, O >::SCoupledElement::node_pt, oomph::FiniteElement::node_pt(), s, SCoupling< M, O >::SCoupledElement::surface_node_position_pt, and oomph::Node::x().

getSolidFeelsParticles()

template<class M , class O >

bool SCoupling< M, O >::getSolidFeelsParticles ( ) const

inline

                                        {
        return solidFeelsParticles_;
    }

setSolidFeelsParticles()

template<class M , class O >

void SCoupling< M, O >::setSolidFeelsParticles ( bool  val )

inline

                                          {
        solidFeelsParticles_ = val;
    }

References calibrate::val.

solveSurfaceCoupling() [1/2]

template<class M , class O >

void SCoupling< M, O >::solveSurfaceCoupling ( const unsigned  max_adapt = 0 )

inline

    {
        // compute nStep
        unsigned nStep = O::getOomphTimeStep()/M::getTimeStep();
        if (nStep==0) {
            // if oomph has a smaller time step than Mercury
            logger(INFO, "Set nStep %, change mercuryTimeStep from % to %",
                   nStep, M::getTimeStep(), O::getOomphTimeStep());
            nStep = 1;
            M::setTimeStep(O::getOomphTimeStep());
        } else {
            logger(INFO, "Set nStep %, change oomphTimeStep from % to %",
                   nStep, O::getOomphTimeStep(), nStep * M::getTimeStep());
            O::setOomphTimeStep(nStep * M::getTimeStep());
        }
 
        // call solve routine
        solveSurfaceCoupling(nStep, max_adapt);
    }

References INFO, and logger.

Referenced by CoupledBeam::CoupledBeam(), and CoupledProblem::CoupledProblem().

solveSurfaceCoupling() [2/2]

template<class M , class O >

void SCoupling< M, O >::solveSurfaceCoupling ( unsigned  nStep, const unsigned  max_adapt  )

inline

    {
        // check whether time steps are set
        logger.assert_always(O::getOomphTimeStep()>0,"Oomph time step not initialised");
        logger.assert_always(M::getTimeStep()>0,"Mercury time step not initialised");
        // check whether isCoupled is set
        logger.assert_always(!coupledBoundaries_.empty(), "isCoupled needs to be set, e.g. via setIsCoupled([](unsigned b) { return b == Boundary::Y_MIN; })");
 
        // first part of the Mercury solve routine, containing the actions before time-stepping
        M::initialiseSolve();
 
        // read Mercury time step
        double mercury_dt = M::getTimeStep();
        logger(INFO, "Mercury time step: %", mercury_dt);
        
        // set oomph time step
        double oomph_dt = nStep * mercury_dt;
        logger(INFO, "Oomph time step: %", oomph_dt);
        
        // set oomph initial conditions
        O::assign_initial_values_impulsive(oomph_dt);
 
        // read oomph-mesh
        logger(INFO, "Set up oomph mesh: % elements (% with traction)",
               O::solid_mesh_pt()?O::solid_mesh_pt()->nelement():0, O::traction_mesh_pt()?O::traction_mesh_pt()->nelement():0);
 
        // get list of bulk elements along the surface-coupled boundaries
        getSCoupledElements();
        
        // create DPM triangle walls from bulk finite elements
        createDPMWallsFromFiniteElems();
        
        // this is the main loop advancing over time
        unsigned nDone = 0; //< last written file number
        while (M::getTime() < M::getTimeMax())
        {
            this->actionsBeforeOomphTimeStep();
            // solve the coupled problem for one time step
            computeOneTimeStepForSCoupling(nStep, max_adapt);
            // write outputs of the oomphProb; this is slaved to the vtk output of Mercury, i.e. an oomph-lib output get written everytime a Mercury vtk file gets written
            if (M::getParticlesWriteVTK() && M::getVtkWriter()->getFileCounter() > nDone)
            {
                O::writeToVTK();
                nDone = M::getVtkWriter()->getFileCounter();
            }
        }
 
        // close output files of mercuryProb
        M::finaliseSolve();
    }

References INFO, and logger.

solveSurfaceCouplingFixedSolid()

template<class M , class O >

void SCoupling< M, O >::solveSurfaceCouplingFixedSolid ( )

inline

    {
        // first part of the Mercury solve routine, containing the actions before time-stepping
        M::initialiseSolve();
        logger.assert_always(!coupledBoundaries_.empty(), "isCoupled needs to be set, e.g. via setIsCoupled([](unsigned b) { return b == Boundary::Y_MIN; })");
        
        // read Mercury time step
        double mercury_dt = M::getTimeStep();
        logger(INFO, "Mercury time step: %", mercury_dt);
        
        // set oomph time step
        logger(INFO, "Solid position fixed");
        
        // set oomph_dt
        this->time_pt()->initialise_dt(0);
        // By default do a non-impulsive start and provide initial conditions
        this->assign_initial_values_impulsive(0);
 
        // read oomph-mesh
        logger(INFO, "Set up oomph mesh: % elements (% with traction)",
               O::solid_mesh_pt()?O::solid_mesh_pt()->nelement():0, O::traction_mesh_pt()?O::traction_mesh_pt()->nelement():0);
        
        // get list of bulk elements along the surface-coupled boundaries
        getSCoupledElements();
        
        // create DPM triangle walls from bulk finite elements
        createDPMWallsFromFiniteElems();
        
        // this is the main loop advancing over time
        unsigned nDone = 0; //< last written file number
        while (M::getTime() < M::getTimeMax())
        {
            this->actionsBeforeOomphTimeStep();
            M::computeOneTimeStep();
            //if (getParticlesWriteVTK() && getVtkWriter()->getFileCounter() > nDone) {
            //    writeToVTK();
            //    nDone = getVtkWriter()->getFileCounter();
            //}
        }
        // close output files of mercuryProb
        M::finaliseSolve();
    }

References INFO, and logger.

solveSurfaceCouplingForgiving()

template<class M , class O >

void SCoupling< M, O >::solveSurfaceCouplingForgiving ( unsigned  nStep, double  timeMaxMin = -constants::inf, const unsigned  max_adapt = 0  )

inline

Solves an unsteady problem, returns successful if timeMaxMin has been reached

                                                                                                                      {
        // solve
        try {
            solveSurfaceCoupling(nStep, max_adapt);
        } catch(OomphLibError& error)  {
            //Store output if newton solver fails
            O::saveSolidMesh();
            M::finaliseSolve();
            double time = O::time_stepper_pt()->time() - nStep * M::getTimeStep();;
            double timeMax = M::getTimeMax();;
            if (time >= timeMaxMin) {
                // take it as successful if a fraction of the time evolution has finished
                logger(INFO,"Newton solver failed at t=% (tMax=%), but code will continue.", time, timeMax);
                exit(0);
            } else {
                logger(ERROR,"Newton solver failed at t=% (tMax=%).", time, timeMax);
            }
        }
    }

References calibrate::error, ERROR, INFO, and logger.

updateDPMWallsFromFiniteElems()

template<class M , class O >

void SCoupling< M, O >::updateDPMWallsFromFiniteElems ( )

inline

    {
        // loop over bulk elements at boundary
        unsigned wallID = 0;
        for (auto sCoupledElement : sCoupledElements_)
        {
            // get members of a SCoupledElement
            Vector<Vector<double*> > position_pt = sCoupledElement.surface_node_position_pt;
            Vector<CoupledSolidNode*> node_pt = sCoupledElement.node_pt;
            
            // reordering vertices of oomph face element (default = {0,1,3,2})
            // \todo can we do this swap in getSCoupledElements, then we don't have to create a local copy position_pt all the time?
            swap(position_pt[2], position_pt[3]);
            
            // get global coordinate at the center
            Vector<double> x(3, 0.0);
            sCoupledElement.bulk_elem_pt->interpolated_x(sCoupledElement.center_local, x);
            Vec3D center {x[0], x[1], x[2]};
            
            // get number of TriangleWalls per oomph face element
            const unsigned nTriangles = position_pt.size();
            unsigned n = 0;
            while (n < nTriangles)
            {
                // get vertices of TriangleWall (multiply vertex position with the length scale of the O)
                std::array<Vec3D, 3> vertex;
                // one vertex at the center
                vertex[0] = center;
                // two vertices from the O<element,TIMESTEPPER>
                vertex[1] = Vec3D(*position_pt[0][0],
                                  *position_pt[0][1],
                                  *position_pt[0][2]);;
                vertex[2] = Vec3D(*position_pt[1][0],
                                  *position_pt[1][1],
                                  *position_pt[1][2]);
                
                // update vertices of triangle facet  (multiply vertex position with the length scale of the O<element,TIMESTEPPER>)
                updateTriangleWall(triangleWalls_[wallID], vertex);
                
                // rotate forward by one element
                rotate(position_pt.begin(), position_pt.begin() + 1, position_pt.end());
                rotate(node_pt.begin(), node_pt.begin() + 1, node_pt.end());
                n++;
                wallID++;
            }
        }
    }

References n, swap(), and plotDoE::x.

updateTractionOnFiniteElems()

template<class M , class O >

void SCoupling< M, O >::updateTractionOnFiniteElems ( )

inline

sets nodal_coupling_residual in each scoupled element

    {
        // if construct mapping with FEM basis functions
        if (!BaseCoupling<M,O>::useCGMapping())
        {
            // tracks the id of the triangle walls (get incremented in computeSCouplingForcesFromTriangles)
            unsigned wallID = 0;
            
            // loop over scoupled elements
            for (auto sCoupledElement : sCoupledElements_)
            {
                // set up memory for nodal coupling force
                Vector<Vector<double> > nodalCouplingForces;
                // returns whether the element is coupled to particles
                bool elemIsCoupled = computeSCouplingForcesFromTriangles(sCoupledElement.bulk_elem_pt,
                                                                         sCoupledElement.surface_node_position_pt.size(),
                                                                         wallID,
                                                                         nodalCouplingForces);
                
                // assign nodal coupling force to the element to be used by element::fill_in_contribution_to_residuals(...)
                sCoupledElement.bulk_elem_pt->set_nodal_coupling_residual(elemIsCoupled, nodalCouplingForces);
            }
            logger(VERBOSE, "Update nodal_coupling_residual");
        }
        else
        {
//            // how many bulk elements in total
//            unsigned n_element = O::solid_mesh_pt()->nelement();
//
//            // loop over the bulk elements
//            for (unsigned e = 0; e < n_element; e++)
//            {
//                // get pointer to the bulk element
//                ELEMENT* elem_pt = dynamic_cast<ELEMENT*>(O::solid_mesh_pt()->element_pt(e));
//
//                // set up memory for nodal coupling force
//                Vector<Vector<double> > nodalCouplingForces;
//                // whether the element is coupled to particles
//                bool elemIsCoupled = computeSCouplingForcesFromCG(elem_pt, nodalCouplingForces);
//
//                // assign nodal coupling force to the element to be used by element::fill_in_contribution_to_residuals(...)
//                elem_pt->set_nodal_coupling_residual(elemIsCoupled, nodalCouplingForces);
//            }
//            // reset the sum of the evaluated CG values
//            for (auto inter : M::interactionHandler)
//            {
//                if (inter->getI()->getName() == "TriangleWall")
//                {
//                    inter->resetTotalCGEval();
//                }
//            }
        }
    }

References logger, and VERBOSE.

updateTriangleWall()

template<class M , class O >

void SCoupling< M, O >::updateTriangleWall ( TriangleWall *&  wall, std::array< Vec3D, 3 >  vertex  )

inline

Update the position of a triangle wall (used in updateDPMWallsFromFiniteElems)

Todo:

why is it using setPrescribedPosition

we need to set velocity

    {
        double time0 = M::getTime();
        double dTime = O::getOomphTimeStep();
        std::array<Vec3D,3> vertex0 = wall->getVertices();
        std::array<Vec3D,3> dVertex = {
            vertex[0] - vertex0[0],
            vertex[1] - vertex0[1],
            vertex[2] - vertex0[2]};
        wall->setPrescribedPosition( [time0, dTime, vertex0, dVertex, wall]( double time ) {
            double f = ( time - time0 ) / dTime;
            std::array<Vec3D, 3> vertex = {
                vertex0[0] + dVertex[0] * f,
                vertex0[1] + dVertex[1] * f,
                vertex0[2] + dVertex[2] * f };
            wall->setVertices( vertex[0], vertex[1], vertex[2] );
            //logger(INFO,"p %",vertex[0]);
            return wall->getPosition();
        } );
        Vec3D velocity = (dVertex[0]+dVertex[1]+dVertex[2])/3./dTime;
        wall->setPrescribedVelocity([velocity] (double time) {
            //logger(INFO,"v %",velocity);
            return velocity;
        });
    }

References f(), BaseInteractable::getPosition(), TriangleWall::getVertices(), BaseInteractable::setPrescribedPosition(), BaseInteractable::setPrescribedVelocity(), TriangleWall::setVertices(), and Jeffery_Solution::velocity().

Member Data Documentation

coupledBoundaries_

template<class M , class O >

std::vector<unsigned> SCoupling< M, O >::coupledBoundaries_

private

Function to determine whether boundary is coupled. Needs to be set before solveSurfaceCoupling is called.

logSurfaceCoupling

template<class M , class O >

bool SCoupling< M, O >::logSurfaceCoupling = true

private

sCoupledElements_

template<class M , class O >

Vector<SCoupledElement> SCoupling< M, O >::sCoupledElements_

private

List of surface-coupled elements.

solidFeelsParticles_

template<class M , class O >

bool SCoupling< M, O >::solidFeelsParticles_ = true

private

Set false for one-way coupling (solid does not feel particles), true for two-way coupling

triangleWalls_

template<class M , class O >

std::vector<TriangleWall*> SCoupling< M, O >::triangleWalls_

private

List of triangle walls used for the surface coupling.

---

The documentation for this class was generated from the following file:

• SCoupling.h