Elmer FEM solver
Elmer is an open source finite element software for multiphysical problems
 All Classes Files Functions Variables Typedefs Macros Groups Pages
particleutils Module Reference
Collaboration diagram for particleutils:
[legend]

Data Types

type  particle_t
 

Public Member Functions

subroutine locateparticleinmeshmarch (ElementIndex, Rinit, Rfin, Init, ParticleStatus, AccurateAtFace, StopFaceIndex, Lambda, Velo, No, ParticleWallKernel)
 
integer, parameter particle_allocated = 1
 
integer, parameter particle_waiting = 2
 
integer, parameter particle_initiated = 3
 
integer, parameter particle_located = 4
 
integer, parameter particle_moving = 5
 
integer, parameter particle_faceboundary = 6
 
integer, parameter particle_partboundary = 7
 
integer, parameter particle_hit = 8
 
integer, parameter particle_ready = 9
 
integer, parameter particle_fixedcoord = 10
 
integer, parameter particle_fixedvelo = 11
 
integer, parameter particle_wallboundary = 12
 
integer, parameter particle_lost = 13
 
integer, parameter particle_ghost = 14
 
type(particle_t), target, save globalparticles
 
subroutine particlestatuscount (Particles)
 
real(kind=dp) function,
dimension(3) 
getparticlecoord (Particles, No)
 
real(kind=dp) function,
dimension(3) 
getparticlevelo (Particles, No)
 
real(kind=dp) function,
dimension(3) 
getparticleforce (Particles, No)
 
subroutine setparticlecoord (Particles, No, Coord)
 
subroutine setparticlevelo (Particles, No, Coord)
 
subroutine setparticleforce (Particles, No, Coord)
 
real(kind=dp) function,
dimension(3) 
getparticleprevcoord (Particles, No)
 
subroutine getparticleuvw (Particles, No, u, v, w)
 
subroutine setparticleuvw (Particles, No, u, v, w)
 
subroutine addparticlecoord (Particles, No, Coord)
 
subroutine addparticlevelo (Particles, No, Coord)
 
subroutine addparticleforce (Particles, No, Coord)
 
integer function getparticlestatus (Particles, No)
 
subroutine setparticlestatus (Particles, No, Status)
 
integer function getparticleelement (Particles, No)
 
subroutine setparticleelement (Particles, No, Index)
 
integer function getparticlenode (Particles, No)
 
subroutine setparticlenode (Particles, No, Index)
 
subroutine markinternalelements (Particles)
 
subroutine setparticlepreliminaries (Particles, dim, TimeOrder)
 
subroutine allocateparticles (Particles, NoParticles)
 
subroutine deletelostparticles (Particles)
 
subroutine eliminateexitingparticles (Particles)
 
subroutine increaseparticles (Particles, NoParticles)
 
subroutine destroyparticles (Particles)
 
subroutine releasewaitingparticles (Particles)
 
integer function changeparticlepartition (Particles)
 
subroutine particlestatistics (Particles, DerOrder)
 
subroutine particleinformation (Particles, ParticleStepsTaken, TimestepsTaken, tottime)
 
real(kind=dp) function characteristicspeed (Particles, No)
 
real(kind=dp) function characteristicelementtime (Particles, No)
 
real(kind=dp) function,
dimension(3) 
randompointinelement (Element, Nodes)
 
subroutine initializeparticles (Particles, InitParticles, AppendParticles)
 
subroutine segmentelementintersection (Mesh, BulkElement, Rinit, Rfin, MinLambda, FaceElement)
 
subroutine segmentelementintersection2 (Mesh, BulkElement, Rinit, Rfin, MinLambda, FaceElement)
 
logical function segmentelementinside (Mesh, BulkElement, Rfin)
 
subroutine locateparticleinmeshoctree (ElementIndex, GlobalCoords, LocalCoords)
 
subroutine locateparticles (Particles, ParticleWallKernel)
 
logical function particleelementinfo (CurrentElement, GlobalCoord, SqrtElementMetric, Basis, dBasisdx)
 
subroutine getvectorfieldinmesh (Var, CurrentElement, Basis, Velo, dBasisdx, GradVelo)
 
subroutine getscalarfieldinmesh (Var, CurrentElement, Basis, Pot, dBasisdx, GradPot)
 
logical function getparticleelementintersection (Particles, BulkElement, Basis, Coord, Radius, BulkElement2, VolumeFraction, AreaFraction)
 
real(kind=dp) function getmaterialpropertyinmesh (PropertyName, BulkElement, Basis, BulkElement2, VolumeFraction)
 
subroutine createneighbourlist (Particles)
 
subroutine createghostparticles (Particles)
 
subroutine destroyghostparticles (Particles)
 
integer function getnextneighbour (Particles, No)
 
subroutine particleparticleinteraction (Particles, dtime, Collision, InteractionKernel)
 
subroutine particleinitializetime (Particles, No)
 
subroutine particleadvancetimestep (Particles, RKstep)
 
subroutine particleboxperiodic (Particles)
 
subroutine particleboxcontact (Particles)
 
real(kind=dp) function getparticletimestep (Particles, InitInterval, tinit)
 
subroutine particlevariablecreate (Particles, Name, DOFs, Output, Secondary, TYPE)
 
type(variable_t) function, pointer particlevariableget (Particles, Name)
 
subroutine particleoutputtable (Particles)
 
subroutine particleoutputgmsh (Particles)
 
subroutine particleoutputvtu (Particles)
 
subroutine particleoutputvti (Particles, GridExtent, GridOrigin, GridDx, GridIndex)
 
subroutine saveparticledata (Model, Solver, dt, TransientSimulation)
 

Member Function/Subroutine Documentation

subroutine particleutils::addparticlecoord ( type(particle_t), pointer  Particles,
integer  No,
real(kind=dp), dimension(3)  Coord 
)

Adds a displacement to the particle coordinates.

subroutine particleutils::addparticleforce ( type(particle_t), pointer  Particles,
integer  No,
real(kind=dp), dimension(3)  Coord 
)

Adds to the force acting on the particle.

Referenced by particlefieldinteraction(), and particleparticleinteraction().

Here is the caller graph for this function:

subroutine particleutils::addparticlevelo ( type(particle_t), pointer  Particles,
integer  No,
real(kind=dp), dimension(3)  Coord 
)

Adds a velocity difference to the particle velocity.

subroutine particleutils::allocateparticles ( type(particle_t), pointer  Particles,
integer  NoParticles 
)

Subroutine allocate particles before launching them.

References messages::info().

Referenced by increaseparticles(), and initializeparticles().

Here is the call graph for this function:

Here is the caller graph for this function:

integer function particleutils::changeparticlepartition ( integer, dimension(:), allocatable, pointer  Particles)

Subroutine for chanching the partition of particles that cross the partition boundary.

References buf, sparitercomm::checkbuffer(), defutils::getmesh(), increaseparticles(), messages::info(), sparitercomm::meshneighbours(), searchelement(), and messages::warn().

Referenced by locateparticles().

Here is the call graph for this function:

Here is the caller graph for this function:

real(kind=dp) function particleutils::characteristicelementtime ( type(particle_t), pointer  Particles,
integer, optional  No 
)

Computes the characterestic time spent in an element Currently computed just for one element as computing the size of element is a timeconsuming operation.

References characteristicspeed(), elementdescription::elementinfo(), elementdescription::elementsize(), defutils::getelementfamily(), defutils::getelementnodes(), defutils::getmesh(), and messages::warn().

Referenced by getparticletimestep().

Here is the call graph for this function:

Here is the caller graph for this function:

real(kind=dp) function particleutils::characteristicspeed ( type(particle_t), pointer  Particles,
integer, optional  No 
)

Computes the characterestic speed for time integration. The speed may be either computed for the whole set or alternatively to just one particle.

References defutils::getlogical(), defutils::getsolverparams(), messages::info(), and parallelutils::parallelreduction().

Referenced by characteristicelementtime(), and getparticletimestep().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::createghostparticles ( integer, dimension(:), allocatable, pointer  Particles)

References buf, sparitercomm::checkbuffer(), defutils::getmesh(), increaseparticles(), messages::info(), sparitercomm::meshneighbours(), and sparitercomm::searchnode().

Referenced by createneighbourlist().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::createneighbourlist ( type(particle_t), pointer  Particles)

This routine creates the nearest neighbours for all nodes. The particle-particle connections may then be found by going through all the nodes of elements.

References createghostparticles(), messages::fatal(), defutils::getelementnodes(), defutils::getelementnofnodes(), and defutils::getmesh().

Referenced by particledynamics().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::deletelostparticles ( type(particle_t), pointer  Particles)

Subroutine deletes lost particles. The reason for losing particles may be that they go to a neighbouring partition.

Referenced by eliminateexitingparticles(), increaseparticles(), and particledynamics().

Here is the caller graph for this function:

subroutine particleutils::destroyghostparticles ( type(particle_t), pointer  Particles)

This assumes that the real particles are followed by ghost particles. which are destroyed before leaving this configuration.

Referenced by particledynamics().

Here is the caller graph for this function:

subroutine particleutils::destroyparticles ( type(particle_t), pointer  Particles)
subroutine particleutils::eliminateexitingparticles ( type(particle_t), pointer  Particles)

Subroutine sets particles that are still sitting on boundary to leave the premises and call for their delition.

References deletelostparticles(), and defutils::getmesh().

Referenced by particledynamics().

Here is the call graph for this function:

Here is the caller graph for this function:

real(kind=dp) function particleutils::getmaterialpropertyinmesh ( character(len=max_name_len)  PropertyName,
type(element_t), pointer  BulkElement,
real(kind=dp), dimension(:)  Basis,
type(element_t), optional, pointer  BulkElement2,
real(kind=dp), optional  VolumeFraction 
)

This subroutine may be used to enquire position dependent material data. Also if the particle is splitted between two elements then this routine can assess the data on the secondary mesh.

References defutils::getmesh(), lists::listgetinteger(), and lists::listgetreal().

Referenced by particlefieldinteraction().

Here is the call graph for this function:

Here is the caller graph for this function:

integer function particleutils::getnextneighbour ( type(particle_t), pointer  Particles,
integer  No 
)

For the first call of given node do the list, thereafter Return the index until the list is finished.

References defutils::getelementnofnodes(), defutils::getmesh(), and messages::info().

Referenced by particleparticleinteraction().

Here is the call graph for this function:

Here is the caller graph for this function:

real(kind=dp) function, dimension(3) particleutils::getparticlecoord ( type(particle_t), pointer  Particles,
integer  No 
)

Returns coordinates of the particle.

Referenced by locateparticles(), particlefieldinteraction(), and particleparticleinteraction().

Here is the caller graph for this function:

integer function particleutils::getparticleelement ( type(particle_t), pointer  Particles,
integer  No 
)

Gets the elements where the particle is located in, or was located last time.

Referenced by locateparticles(), particlefieldinteraction(), setadvectedfield(), and setparticlevelocities().

Here is the caller graph for this function:

logical function particleutils::getparticleelementintersection ( type(particle_t), pointer  Particles,
type(element_t), pointer  BulkElement,
real(kind=dp), dimension(:)  Basis,
real(kind=dp), dimension(3)  Coord,
real(kind=dp)  Radius,
type(element_t), pointer  BulkElement2,
real(kind=dp)  VolumeFraction,
real(kind=dp), optional  AreaFraction 
)

The routine returns the possible intersection of a secondary element with a different material property and the circle / sphere. For example, the buoyancy at the interface will depend on the weighted sum of the densities of the two materials.

References defutils::getelementnodes(), defutils::getmesh(), lists::listgetinteger(), normal(), and elementdescription::pointfacedistance().

Referenced by particlefieldinteraction().

Here is the call graph for this function:

Here is the caller graph for this function:

real(kind=dp) function, dimension(3) particleutils::getparticleforce ( type(particle_t), pointer  Particles,
integer  No 
)

Returns force acting on the particle.

integer function particleutils::getparticlenode ( type(particle_t), pointer  Particles,
integer  No 
)

Gets the closest node related to the particle.

real(kind=dp) function, dimension(3) particleutils::getparticleprevcoord ( type(particle_t), pointer  Particles,
integer  No 
)

Gets the previous particle coordinate.

Referenced by locateparticles().

Here is the caller graph for this function:

integer function particleutils::getparticlestatus ( type(particle_t), pointer  Particles,
integer  No 
)

Gets the status of the particle.

Referenced by particlefieldinteraction(), setadvectedfield(), and setparticlevelocities().

Here is the caller graph for this function:

real(kind=dp) function particleutils::getparticletimestep ( type(particle_t), pointer  Particles,
logical  InitInterval,
real(kind=dp), optional  tinit 
)

Set a the timestep for the particles. Depending on the definitions the timestep may be the same for all particles, or may be defined independently for each particle.

References characteristicelementtime(), characteristicspeed(), messages::fatal(), defutils::getcreal(), defutils::getinteger(), defutils::getlogical(), defutils::getsolverparams(), defutils::gettimestepsize(), particlevariablecreate(), and particlevariableget().

Referenced by particleadvector(), and particledynamics().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::getparticleuvw ( type(particle_t), pointer  Particles,
integer  No,
real(kind=dp)  u,
real(kind=dp)  v,
real(kind=dp), optional  w 
)

Gets the local coordinates of the element of the given particle.

real(kind=dp) function, dimension(3) particleutils::getparticlevelo ( type(particle_t), pointer  Particles,
integer  No 
)

Returns velocity of the particle.

Referenced by locateparticles(), particlefieldinteraction(), and particleparticleinteraction().

Here is the caller graph for this function:

subroutine particleutils::getscalarfieldinmesh ( type(variable_t), pointer  Var,
type(element_t)  CurrentElement,
real(kind=dp), dimension(:)  Basis,
real(kind=dp)  Pot,
real(kind=dp), dimension(:,:), optional  dBasisdx,
real(kind=dp), dimension(:), optional  GradPot 
)

The routine returns a potential and its gradient.

References defutils::getmesh().

Referenced by particlefieldinteraction(), and setadvectedfield().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::getvectorfieldinmesh ( type(variable_t), pointer  Var,
type(element_t)  CurrentElement,
real(kind=dp), dimension(:)  Basis,
real(kind=dp), dimension(:)  Velo,
real(kind=dp), dimension(:,:), optional  dBasisdx,
real(kind=dp), dimension(:,:), optional  GradVelo 
)

The routine returns velocity and optionally a gradient of velocity. These kind of functions are needed repeated and therefore to reduced the size of individual solvers it has been hard coded here.

References defutils::getlogical(), defutils::getmesh(), and defutils::getsolverparams().

Referenced by particlefieldinteraction(), particlewallcontact(), setadvectedfield(), and setparticlevelocities().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::increaseparticles ( type(particle_t), pointer  Particles,
integer  NoParticles 
)

Increase particle array size by given amount.

References allocateparticles(), and deletelostparticles().

Referenced by changeparticlepartition(), and createghostparticles().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::initializeparticles ( type(particle_t), pointer  Particles,
integer, optional  InitParticles,
logical, optional  AppendParticles 
)
subroutine particleutils::locateparticleinmeshmarch ( integer  ElementIndex,
real(kind=dp), dimension(3)  Rinit,
real(kind=dp), dimension(3)  Rfin,
logical  Init,
integer  ParticleStatus,
logical  AccurateAtFace,
integer  StopFaceIndex,
real(kind=dp)  Lambda,
real(kind=dp), dimension(3)  Velo,
integer, optional  No,
optional  ParticleWallKernel 
)

Locate the particle using controlled marching from element to element. The crossing point between given trajectory and all face elements is computed. The one that is passed at first is associated to the next bulk element.

References defutils::getelementnodes(), defutils::getelementnofnodes(), defutils::getmesh(), defutils::getsolverparams(), lists::listgetconstreal(), lists::listgetinteger(), lists::listgetlogical(), normal(), elementdescription::normalvector(), interpolation::pointinelement(), segmentelementinside(), segmentelementintersection(), segmentelementintersection2(), and messages::warn().

Referenced by locateparticles().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::locateparticleinmeshoctree ( integer  ElementIndex,
real(kind=dp), dimension(3)  GlobalCoords,
real(kind=dp), dimension(3), optional  LocalCoords 
)

Find the particle in the mesh using actree based search. This could be preferred in the initial finding of the correct elements. The major downside of the method is that there is no controlled face detection needed for wall interaction, for example.

References interpolation::buildquadranttree(), interpolation::findleafelements(), defutils::getelementnodes(), defutils::getelementnofnodes(), defutils::getmesh(), interpolation::pointinelement(), and messages::warn().

Referenced by asciipointstomesh(), and initializeparticles().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::locateparticles ( type(particle_t), pointer  Particles,
optional  ParticleWallKernel 
)

Locate the particles in their new positions in the mesh.

References changeparticlepartition(), defutils::getmesh(), getparticlecoord(), getparticleelement(), getparticleprevcoord(), getparticlevelo(), defutils::getsolverparams(), lists::listgetlogical(), locateparticleinmeshmarch(), setparticlecoord(), setparticlevelo(), and messages::warn().

Referenced by particleadvector(), and particledynamics().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::markinternalelements ( type(particle_t), pointer  Particles)

The subroutine marks the elements which are not on the boundary, either internal or external one. This information may be used to speed up different loops where particle-boundary interaction is needed.

References messages::fatal(), defutils::getmesh(), and lists::listgetinteger().

Referenced by setparticlepreliminaries().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::particleadvancetimestep ( type(particle_t), pointer  Particles,
integer, optional  RKstep 
)

Advance the particles with a time step. The timestep may also be an intermediate Runge-Kutta step.

References messages::fatal(), defutils::getcreal(), defutils::getsolverparams(), lists::listgetconstreal(), and particlevariableget().

Referenced by particleadvector(), and particledynamics().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::particleboxcontact ( type(particle_t), pointer  Particles)

Checks the boundaries for rectangular and hexahedral shapes and enforces elastic reflection. This is alternative and computationally more economic way and is ideal for testing purposes, at least.

References messages::fatal(), defutils::getcreal(), defutils::getmesh(), defutils::getsolverparams(), lists::listgetintegerarray(), and lists::listgetlogical().

Referenced by particledynamics().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::particleboxperiodic ( type(particle_t), pointer  Particles)

Checks the boundaries for rectangular and hexahedral shapes and enforces periodic BCs. Currently the only supported way for setting periodic BCs.

References defutils::getmesh(), defutils::getsolverparams(), lists::listgetintegerarray(), lists::listgetlogical(), and lists::listgetstring().

Referenced by particledynamics().

Here is the call graph for this function:

Here is the caller graph for this function:

logical function particleutils::particleelementinfo ( type(element_t), pointer  CurrentElement,
real(kind=dp), dimension(:)  GlobalCoord,
real(kind=dp)  SqrtElementMetric,
real(kind=dp), dimension(:)  Basis,
real(kind=dp), dimension(:,:), optional  dBasisdx 
)

Given the element & global coordinates returns the local coordinates. The idea of this routine is to transparently block the local coordinate search from the user by directly giving the basis function values related to a global coordinate. Sloppy tolerances are used since we should have already located the element.

References elementdescription::elementinfo(), defutils::getelementnodes(), interpolation::pointinelement(), and messages::warn().

Referenced by particlefieldinteraction(), particlewallcontact(), setadvectedfield(), and setparticlevelocities().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::particleinformation ( type(particle_t), pointer  Particles,
integer  ParticleStepsTaken,
integer  TimestepsTaken,
real(kind=dp)  tottime 
)

Echos some information on particle amounts to standard output.

References messages::info(), parallelutils::parallelreduction(), and particlestatuscount().

Referenced by particleadvector(), and particledynamics().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::particleinitializetime ( type(particle_t), pointer  Particles,
integer, optional  No 
)

Initialize the time for next time integration step.

Referenced by particledynamics().

Here is the caller graph for this function:

subroutine particleutils::particleoutputgmsh ( type(particle_t), pointer  Particles)

Write particles to an external file in Gmsh format as vector points.

References defutils::getmesh(), defutils::getsolverparams(), messages::info(), lists::listgetintegerarray(), lists::listgetstring(), and lists::variableget().

Referenced by saveparticledata().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::particleoutputtable ( type(particle_t), pointer  Particles)

Write particles to an external file in simple ascii format matrix.

References closeparticlefile(), messages::fatal(), defutils::getmesh(), defutils::getsolverparams(), lists::listgetinteger(), lists::listgetlogical(), lists::listgetstring(), openparticlefile(), particlevariableget(), lists::variableget(), writeparticlefilenames(), and writeparticleline().

Referenced by savegriddata(), and saveparticledata().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::particleoutputvti ( type(particle_t), pointer  Particles,
integer, dimension(6)  GridExtent,
real(kind=dp), dimension(3)  GridOrigin,
real(kind=dp), dimension(3)  GridDx,
integer, dimension(:,:,:), pointer  GridIndex 
)

Writes data out in XML VTK ImageData format (VTI) which assumes a uniform grid where the position of each point is defined by the origin and the grid density. Also binary output and single precision therein is supported.

References defutils::getlogical(), defutils::getmesh(), defutils::getsolverparams(), messages::info(), lists::listgetinteger(), lists::listgetlogical(), lists::listgetstring(), messages::warn(), writepvtifile(), and writevtifile().

Referenced by savegriddata().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::particleoutputvtu ( type(particle_t), pointer  Particles)

Saves particles in unstructured XML VTK format (VTU) to an external file.

References defutils::getmesh(), defutils::getsolverparams(), messages::info(), lists::listgetcreal(), lists::listgetinteger(), lists::listgetlogical(), lists::listgetstring(), writepvtufile(), and writevtufile().

Referenced by savegriddata(), and saveparticledata().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::particleparticleinteraction ( type(particle_t), pointer  Particles,
real(kind=dp)  dtime,
logical  Collision,
  InteractionKernel 
)

Computes interaction between particles given an interaction kernel that that is the pointer to the function that computes the effect of the interaction in terms of forces or new coordinate values (for collision models).

References addparticleforce(), getnextneighbour(), getparticlecoord(), getparticlevelo(), and setparticlecoord().

Referenced by particledynamics().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::particlestatistics ( type(particle_t), pointer  Particles,
integer  DerOrder 
)

Subroutine computes the means of coordinates / velocities / force.

References messages::fatal(), parallelutils::parallelreduction(), and messages::warn().

Referenced by particledynamics().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::particlestatuscount ( type(particle_t), pointer  Particles)

Counts particles in different categories.

References messages::info().

Referenced by particleinformation().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::particlevariablecreate ( type(particle_t), pointer  Particles,
character(len=*)  Name,
integer, optional  DOFs,
logical, optional  Output,
logical, optional  Secondary,
integer, optional  TYPE 
)

Creates a variable related to the particles. The normal variabletype without the permutation vector is used.

References messages::info(), solver(), lists::variableadd(), and messages::warn().

Referenced by getparticletimestep(), particleadvector(), particledynamics(), and particledynamicsstuff::particlewallproc().

Here is the call graph for this function:

Here is the caller graph for this function:

type(variable_t) function, pointer particleutils::particlevariableget ( type(particle_t), pointer  Particles,
character(len=*)  Name 
)

Given a variable name, get a handle to it.

References lists::variableget().

Referenced by getparticletimestep(), particleadvancetimestep(), particleadvector(), particlefieldinteraction(), particleoutputtable(), particledynamicsstuff::particlewallproc(), setadvectedfield(), setfixedparticles(), setparticlevelocities(), and writevtufile().

Here is the call graph for this function:

Here is the caller graph for this function:

real(kind=dp) function, dimension(3) particleutils::randompointinelement ( type(element_t)  Element,
type(nodes_t)  Nodes 
)

Finds a random point that is guaranteed within the element.

References elementdescription::elementinfo(), generalutils::evenrandom(), and messages::fatal().

Referenced by initializeparticles().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::releasewaitingparticles ( type(particle_t), pointer  Particles)

Subroutine for releaseing initiated but waiting particles.

References defutils::getcreal(), defutils::getinteger(), and defutils::getsolverparams().

Referenced by particleadvector(), and particledynamics().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::saveparticledata ( type(model_t)  Model,
type(solver_t), target  Solver,
real(kind=dp)  dt,
logical  TransientSimulation 
)

Write particles to an external file in various formats.

References defutils::getsolverparams(), lists::listaddstring(), lists::listcheckpresent(), lists::listgetlogical(), lists::listgetstring(), particleoutputgmsh(), particleoutputtable(), particleoutputvtu(), and messages::warn().

Referenced by particledynamics(), and particleoutputsolver().

Here is the call graph for this function:

Here is the caller graph for this function:

logical function particleutils::segmentelementinside ( type(mesh_t), pointer  Mesh,
type(element_t), pointer  BulkElement,
real(kind=dp), dimension(3)  Rfin 
)

This subroutine tests whether a particle is within element using the consistant strategy with the above algorithm.

References defutils::getelementnodes(), defutils::getelementnofnodes(), and elementdescription::linefaceintersection2().

Referenced by locateparticleinmeshmarch().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::segmentelementintersection ( type(mesh_t), pointer  Mesh,
type(element_t), pointer  BulkElement,
real(kind=dp), dimension(3)  Rinit,
real(kind=dp), dimension(3)  Rfin,
real(kind=dp)  MinLambda,
type(element_t), pointer  FaceElement 
)

This subroutine finds the possible intersection between elementfaces and a line segment defined by two coordinates.

References defutils::getelementnodes(), elementdescription::linefaceintersection(), generalutils::sortr(), and messages::warn().

Referenced by locateparticleinmeshmarch().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::segmentelementintersection2 ( type(mesh_t), pointer  Mesh,
type(element_t), pointer  BulkElement,
real(kind=dp), dimension(3)  Rinit,
real(kind=dp), dimension(3)  Rfin,
real(kind=dp)  MinLambda,
type(element_t), pointer  FaceElement 
)

This subroutine finds the possible intersection between elementfaces and a line segment defined by two coordinates.

References defutils::getelementnodes(), and elementdescription::linefaceintersection2().

Referenced by locateparticleinmeshmarch().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::setparticlecoord ( type(particle_t), pointer  Particles,
integer  No,
real(kind=dp), dimension(3)  Coord 
)

Sets the particle coordinates.

Referenced by locateparticles(), and particleparticleinteraction().

Here is the caller graph for this function:

subroutine particleutils::setparticleelement ( type(particle_t), pointer  Particles,
integer  No,
integer  Index 
)

Sets the element where the particle is located.

subroutine particleutils::setparticleforce ( type(particle_t), pointer  Particles,
integer  No,
real(kind=dp), dimension(3)  Coord 
)

Sets the particle force.

subroutine particleutils::setparticlenode ( type(particle_t), pointer  Particles,
integer  No,
integer  Index 
)

Sets the closest node related to the particle.

subroutine particleutils::setparticlepreliminaries ( type(particle_t), pointer  Particles,
integer, optional  dim,
integer, optional  TimeOrder 
)

Subroutine sets up some preliminary information needed for the particler tracker: timeorder, space dimension, bounding box, and mesh edges/faces.

References messages::fatal(), meshutils::findmeshedges(), defutils::getmesh(), markinternalelements(), sparitercomm::sparedgenumbering(), and sparitercomm::sparfacenumbering().

Referenced by particleadvector(), and particledynamics().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine particleutils::setparticlestatus ( type(particle_t), pointer  Particles,
integer  No,
integer  Status 
)

Sets the status of the particle.

Referenced by particlefieldinteraction().

Here is the caller graph for this function:

subroutine particleutils::setparticleuvw ( type(particle_t), pointer  Particles,
integer  No,
real(kind=dp)  u,
real(kind=dp)  v,
real(kind=dp), optional  w 
)

Sets the local coordinates of the element of the given particle.

subroutine particleutils::setparticlevelo ( type(particle_t), pointer  Particles,
integer  No,
real(kind=dp), dimension(3)  Coord 
)

Sets the particle velocity.

Referenced by locateparticles().

Here is the caller graph for this function:

Member Data Documentation

type(particle_t), target, save particleutils::globalparticles
integer, parameter particleutils::particle_allocated = 1
integer, parameter particleutils::particle_faceboundary = 6
integer, parameter particleutils::particle_fixedcoord = 10
integer, parameter particleutils::particle_fixedvelo = 11
integer, parameter particleutils::particle_ghost = 14
integer, parameter particleutils::particle_hit = 8
integer, parameter particleutils::particle_initiated = 3
integer, parameter particleutils::particle_located = 4
integer, parameter particleutils::particle_lost = 13
integer, parameter particleutils::particle_moving = 5
integer, parameter particleutils::particle_partboundary = 7
integer, parameter particleutils::particle_ready = 9
integer, parameter particleutils::particle_waiting = 2
integer, parameter particleutils::particle_wallboundary = 12

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