migflow.fluid

Model for Immersed Granular Flow: Fluid user interface

Contact: jonathan.lambrechts@uclouvain.be Webpage: www.migflow.be

MigFlow computes immersed granular flows using an unresolved FEM-DEM model. The continuous phase (mixture) is computed by solving averaged Navier-Stokes equations on unstructured meshes with the continuous finite element method.

class migflow.fluid.FluidProblem(dim, g, mu, rho, sigma=0, coeff_stab=0.01, volume_drag=0.0, quadratic_drag=0.0, solver=None, solver_options='', drag_in_stab=1, drag_coefficient_factor=1, temporal=True, advection=True, usolid=False, p2p1=False, full_implicit=False, model_b=False)

Bases: object

Creates the numerical representation of the fluid.

Builds the fluid structure.

Parameters:
  • dim (int) – Dimension of the problem
  • g (ndarray) – Gravity vector
  • mu (list) – List of fluid phases dynamic viscosities (this should be a vector whom dimension specify if there is one or two fluid)
  • rho (list) – List of fluid phases densities (this should be a vector whom dimension specify if there is one or two fluid)
  • sigma (float) – Surface tension (only when there are two fluids)
  • coeff_stab (float) – Optional argument used as coefficient for extra diffusivity added to stabilise the advection of concentration (only for two fluid problem)
  • volume_drag (float) – Enables volumetric drag : -volume_drag * u[i]
  • quadratic_drag (float) – Enables quadratic drag : -quadratic_drag*||uold|| * u[i]
  • solver (Optional[str]) – possible solver are “mumps”, “petsc”, “petsc4py”, “scipy”
  • solver_options (str) – Optional argument to specify solver option
  • drag_in_stab (int) – States if the drag force is in the stabilisation term
  • drag_coefficient_factor (int) – Factor multiplicating the drag coefficient
  • temporal (bool) – Enables d/dt (i.e. False = steady)
  • advection (bool) – Enables advective terms (i.e. False = Stokes, True = Navier-Stokes
  • usolid (bool) – if False, div.u_solid is evaluated = (porosity-old_porosity)/dt
  • p2p1 (bool) – if True, use P2 element for velocity fields and P1 element for pressure field, else a stabilised P1-P1 formulation is used
  • full_implicit (bool) – If False, the linearised averaged Navier-Stokes equation are considered else a fully implicit formulation is used
  • model_b (bool) – If True, the buoyancy term is assumed to be expressed as an additional drag
Raises:
  • ValueError – If the dimension differs from 2 or 3
  • NameError – If C builder cannot be found
__del__()

Deletes the fluid structure.

adapt_mesh(path_to_mesh='adapt/mesh', lcmax=None, lcmin=None, n_el=None, old_n_particles=0, old_particle_position=None, old_particle_volume=None)

Adapts the mesh.

Parameters:
  • old_n_particles (int) – Number of particles at the previous time step
  • old_particle_position (Optional[ndarray]) – Position of the particles at the previous time step
  • old_particle_volume (Optional[ndarray]) – Volume of the particles at the previous time step
add_constraint(constraint_index, constraint_value)
advance_concentration(dt)

Solves the advection equation for the concentration using the current velocity field.

Parameters:dt (float) – Computation time step
boundary_force()

Give access to force exerted by the fluid on the boundaries.

Returns:Forces exerted by the fluid on the boundaries.
Return type:np.ndarray
bulk_force()

Gives access to the volume force at fluid nodes. :returns: Bulk force :rtype: np.ndarray

compute_cfl()

Computes the CFL number divided by the time step

compute_node_force(dt)

Computes the forces to apply on each grain resulting from the fluid motion.

Parameters:dt (float) – Computation time step
compute_node_torque(dt)

Computes the angular drags to apply on each grain resulting from the fluid motion. Only in 2D :type dt: float :param dt: Computation time step

concentration_dg()

Returns the concentration at discontinuous nodes.

Returns:concentration field
Return type:np.ndarray
concentration_dg_grad()

Returns the gradient of the concentration field at discontinuous nodes.

Returns:concentration gradient
Return type:np.ndarray
coordinates()

Gives access to the coordinate of the mesh nodes.

Returns:mesh nodal coordinates
Return type:np.ndarray
coordinates_fields()

Returns the coordinates of each degree of freedom of the solution.

Returns:fields coordinates
Return type:np.ndarray
curvature()

Returns the curvature at each element.

Returns:curvature
Return type:np.ndarray
dimension()

Returns the dimension of the problem

Returns:dimension
Return type:int
elements()

Gives read-only access to the elements of the mesh.

Returns:mesh elements
Return type:np.ndarray
field_index(ifield)

Gives access to all the index of a field into the solution vector.

Parameters:ifield (int) – field index (0,..,dimension-1) is a velocity, dimension gives the pressure.
Returns:index array.
Return type:np.ndarray
fields_gradient()

Returns an continuous gradient estimation of each field on a P1 mesh.

Returns:Fields_gradient
Return type:np.ndarray
full_implicit_euler(**kwargs)
get_default_export(dtype=<class 'numpy.float64'>)

Returns the default fields to write in your output as a dictionnary {field : (values, discretisation_element)} :param dtype: numpy representation. Defaults to np.float64.

Returns:default fields
Return type:dict
get_mapping(etype)

Get mapping associated to an element.

Parameters:etype (str) – element type (“P1”, “P1DG”, “P2”, “P2DG”)
Returns:mapping associated to the element type.
Return type:np.ndarray
get_p1_element()

Returns the P1 element used -> “triangle_p1” or “tetrahedron_p1”

Returns:“triangle_p1”
Return type:str
get_pressure_element()

Returns the discretisation element of the pressure field :returns: pressure element :rtype: str

get_velocity_element()

Returns the discretisation element of the velocity field :returns: velocity element :rtype: str

global_map()

Gives access to the map of all the degrees of freedom.

Returns:global mapping
Return type:np.ndarray
implicit_euler(**kwargs)
interpolate(solution=None, velocity=None, velocity_x=None, velocity_y=None, velocity_z=None, pressure=None)

Imposes the solution (or the field) to the prescribed value. If a callback is given, the solution (or field) is estimated at its local position. The solution (or field) is strongly imposed (not weakly).

Parameters:
  • solution (Union[ndarray, <built-in function callable>, None]) – imposed solution vector
  • velocity (Union[ndarray, <built-in function callable>, None]) – imposed velocity vector.
  • velocity_x (Union[ndarray, <built-in function callable>, None]) – imposed horizontal velocity component.
  • velocity_y (Union[ndarray, <built-in function callable>, None]) – imposed horizontal velocity component.
  • velocity_z (Union[ndarray, <built-in function callable>, None]) – imposed horizontal velocity component.
  • pressure (Union[ndarray, <built-in function callable>, None]) – imposed pressure field.
load_msh(filename)

Sets the domain geometry for the fluid computation.

Parameters:filename (str) – Name of the msh file. If None the current gmsh model is loaded without the need to a msh file.
local_size()

Returns the number of degree of freedom by element.

Returns:number of degree of freedom by element
Return type:int
mesh_boundaries()

Returns the mesh boundaries information as a dictionnary : {boundary_number : edges} :returns: Mesh information : {boundary_number : edges} :rtype: dict

mesh_velocity()

Give access to the mesh velocity value at the mesh nodes.

Returns:mesh velocities at the mesh nodes
Return type:np.ndarray
move_particles(position, velocity, contact)

Set location of the grains in the mesh and compute the porosity in each cell.

Parameters:
  • position (ndarray) – List of particles centre positions
  • velocity (ndarray) – List of particles velocity
  • contact (ndarray) – List of particles contact resultant forces
n_elements()

Returns the number of mesh elements.

Returns:number of elements
Return type:int
n_fields()

Returns the number of fluid fields.

Returns:number of fields
Return type:int
n_fluids()

Returns the number of fluids (only one or two).

Returns:Number of fluids
Return type:int
n_nodes()

Returns the number of mesh nodes.

Returns:number of nodes
Return type:int
node_volume()

Returns the volume associated with each node. :returns: node volume :rtype: np.ndarray

old_porosity()

Returns the old porosity.

Returns:old porosity
Return type:np.ndarray
parent_nodes()

Gives access to the parent nodes of each node.

Returns:The parent nodes mapping
Return type:np.ndarray
particle_element_id()

Returns the id of the mesh cell in which particles are located.

Returns:particle element id
Return type:np.ndarray
particle_position()

Gives access to the particles position. :returns: particles position :rtype: np.ndarray

particle_uvw()

Returns the coordinates of the particles inside their element :returns: parametric coordinates of the particles :rtype: np.ndarray

particle_velocity()

Gives access to the particles velocity. :returns: particles velocity :rtype: np.ndarray

particle_volume()

Gives access to the particles volume. :returns: particles volume :rtype: np.ndarray

porosity()

Returns the porosity at independant nodes.

Returns:volume fluid fraction
Return type:np.ndarray
pressure()

Reads the pressure solution.

Returns:pressure
Return type:np.ndarray
pressure_index()

Returns the index of the pressure field into the solution vector.

Returns:pressure index array
Return type:np.ndarray
read_mig(odir, t)

Reads output file to reload computation.

Parameters:
  • odir (str) – Directory in which to read the file
  • t (float) – Time at which to read the file
read_vtk(dirname, i)

(DEPRECATED) Reads output file to reload computation.

Parameters:filename – Name of the file to read
set_concentration_cg(concentration)

Sets the concentration at nodes. :type concentration: ndarray :param concentration: concentration field

set_coordinates(x)

Sets the coordinates of the mesh nodes

Parameters:x (ndarray) – new nodes positions
set_mean_pressure(mean_pressure)

Enables the nodal mean pressure constraint and imposes its value.

Parameters:mean_pressure (float) – add a constraint to impose sum_i vol_i p_i = mean_pressure*volume
set_mesh(nodes, elements, boundaries)

Sets the domain geometry for the fluid computation.

Parameters:
  • nodes (ndarray) – (n_node, dimension) numpy float array of the nodes coordinates
  • elements (ndarray) – (n_elements, dimension+1) numpy int array of the elements
  • boundaries (list) – [(name, (n_elements, dimension))] list of pair of name (string)
  • numpy int array with the corresponding boundary elements (and) –
set_open_boundary(tag, velocity=None, pressure=None, concentration=None, viscous_flag=1)

Sets the weak open boundaries (=normal fluxes) for the fluid problem.

Parameters:
  • tag (Union[str, List[str]]) – Physical tag (or list of tags), set in the mesh, of the open boundaries
  • velocity (Union[float, <built-in function callable>, None]) – The velocity value if imposed (callback or number)
  • pressure (Union[float, <built-in function callable>, None]) – The pressure value if imposed (callback or number)
  • concentration (Union[float, <built-in function callable>, None]) – Concentration outside the boundary
  • viscous_flag (bool) – Flag to compute the viscous term at the boundary
set_particles(delassus, volume, position, velocity, contact)

Set location of the grains in the mesh and compute the porosity in each cell.

Parameters:
  • delassus (ndarray) – List of particles delassus operators
  • volume (ndarray) – List of particles volume
  • position (ndarray) – List of particles centre positions
  • velocity (ndarray) – List of particles velocity
  • contact (ndarray) – List of particles contact resultant forces
set_strong_boundary(tag, pressure=None, velocity=None, velocity_x=None, velocity_y=None, velocity_z=None, solution=None)

Sets the strong boundaries (=constrains) for the fluid problem.

Parameters:
  • tag (str) – Physical tag (set in the mesh) of the boundary on which the constraint is added
  • pressure (Union[float, <built-in function callable>, None]) – value or callback assigned to the pressure field
  • velocity (Union[float, <built-in function callable>, None]) – value or callback assigned to the velocity field
  • velocity_x (Union[float, <built-in function callable>, None]) – value or callback assigned to the velocity_x field
  • velocity_y (Union[float, <built-in function callable>, None]) – value or callback assigned to the velocity_y field
  • velocity_z (Union[float, <built-in function callable>, None]) – value or callback assigned to the velocity_z field
  • solution (Union[float, <built-in function callable>, None]) – value or callback assigned to the solution field
set_symmetry_boundary(tag, pressure=None)

Sets the weak symmetry boundaries (=normal fluxes) for the fluid problem.

Parameters:
  • tag (Union[str, List[str]]) – Physical tag (or list of tags), set in the mesh file, of the symmetry boundaries
  • pressure (Union[float, <built-in function callable>, None]) – The pressure value if imposed (callback or number)
set_wall_boundary(tag, pressure=None, velocity=None, viscous_flag=None)

Sets the weak wall boundaries (=normal fluxes) for the fluid problem.

Parameters:
  • tag (Union[str, List[str]]) – Physical tag (or list of tags), set in the mesh, of the wall boundaries
  • pressure (Union[float, <built-in function callable>, None]) – The pressure value if imposed (callback or number)
  • velocity (Union[float, <built-in function callable>, None]) – The velocity value if imposed (callback or number)
  • viscous_flag (Optional[bool]) – Flag to compute the viscous term at the boundary
set_weak_boundary(tag, pressure=None, velocity=None, concentration=None, viscous_flag=None)

Sets the weak boundaries (=normal fluxes) for the fluid problem.

Parameters:
  • tag (Union[str, List[str]]) – Physical tag (or list of tags), set in the mesh, of the weak boundaries
  • pressure (Union[float, <built-in function callable>, None]) – The pressure value if imposed (callback or number)
  • velocity (Union[float, <built-in function callable>, None]) – The velocity value if imposed (callback or number)
  • concentration (Union[float, <built-in function callable>, None]) – Concentration outside the boundary
  • viscous_flag (Optional[bool]) – Flag to compute the viscous term at the boundary
set_weak_boundary_nodal_values(tag, values)

Set the nodal values to use for the weak boundary instead of using a callback or a prescribed value. Only available for p1p1 formulation.

Parameters:
  • tag (str) – Physical tag (or list of tags), set in the mesh
  • values (ndarray) – nodal values imposed on the boundary
Raises :
ValueError: if the tag is not found
solution()

Gives access to the fluid field value at each degree of freedom as a flat array.

Returns:fluid solution
Return type:np.ndarray
solution_at_coordinates(x)

Returns the solution vector interpolated at the given coordinates.

Parameters:x (ndarray) – coordinates
Returns:solution at given coordinates
Return type:np.ndarray
u_solid()

Returns the solid velocity as a p1-continuous field.

Returns:solid velocity
Return type:np.ndarray
velocity()

Reads the velocity solution.

Returns:velocity
Return type:np.ndarray
velocity_index()

Returns the index of the velocities into the solution vector as a bidimensionnal array.

Returns:velocities index array
Return type:np.ndarray
volume_flux(bnd_tag)

Computes the integral of the (outgoing) normal velocity through boundary with tag bnd_tag.

Returns:volume flux
Return type:float
write_mig(output_dir, t, fields=None, mesh_dtype=<class 'numpy.float64'>)

Writes the output files for post-visualization. Metadata are stored into a file called “fluid.migf” while the data are stored into the “fluid” directory. :type output_dir: str :param output_dir: Output directory :type t: float :param t: Computational time :param extra_fields: extra field as a dictionnary {name : values at p1 degree of freedom}

write_vtk(output_dir, it, t)

(DEPRECATED) Writes output file for post-visualization. :param output_dir: Output directory :param it: Computation iteration :param t: Computation time

migflow.fluid.removeprefix(str, prefix)
Return type:str