medypt.utils ============ .. py:module:: medypt.utils .. autoapi-nested-parse:: Convenient tools for the package. Attributes ---------- .. autoapisummary:: medypt.utils.SQRT_PI Classes ------- .. autoapisummary:: medypt.utils.L2RError medypt.utils.Projector medypt.utils.TMat Functions --------- .. autoapisummary:: medypt.utils.create_mesh medypt.utils.f1_2 medypt.utils.d_density medypt.utils.ufl_mat2voigt4strain medypt.utils.young_poisson2stiffness medypt.utils.ufl_tr_voigt medypt.utils.gen_therm_noise Module Contents --------------- .. py:data:: SQRT_PI :value: 1.77245385091 Square root of pi. .. py:function:: create_mesh(comm: mpi4py.MPI.Comm, model: gmsh.model | str, out_file: str | None = None, **kwargs: Any) -> dolfinx.io.gmsh.MeshData Create a DOLFINx mesh from a Gmsh model and output to file. :param comm: MPI communicator top create the mesh on. :type comm: MPI.Comm :param model: Gmsh model or name of a file containing the mesh (``.msh`` or ``.xdmf``). Reading from ``.xdmf`` file is the most efficient and is recommended for large mesh. :type model: gmsh.model | str :param out_file: ``.xdmf`` file name for writing. Default to ``None`` for not writing to ``.xdmf`` file. :type out_file: str | None :param kwargs: Additional keyword arguments: * mesh_dim (int): Geometric dimension of mesh of the gmsh model or ``.msh`` file. Required if ``model`` is a :py:class:`gmsh.model` or a ``.msh`` file. * rank (int): Rank of the MPI process used for generating from gmsh model or reading from ``.msh`` files. Required if ``model`` is a :py:class:`gmsh.model` or a ``.msh`` file. * mesh_name (str): Name (identifier) of the mesh to read from the input ``.xdmf`` file and to add to the output file. Required if ``model`` is a ``.xdmf`` file. * mode (str): Mode for writing mesh to ``.xdmf`` file. ``'w'`` (write) or ``'a'`` (append). Required if ``out_file`` is not ``None``. :type kwargs: dict[str, Any] :returns: The created mesh data. :rtype: gmshio.MeshData .. attention:: This function is copied from https://docs.fenicsproject.org/dolfinx/main/python/demos/demo_gmsh.html and modified. .. py:function:: f1_2(x: Any, exp: collections.abc.Callable = ufl.exp) -> Any Approximate analytic expression for Fermi-Dirac integral of order 1/2. :param x: Input value. :type x: Any :param exp: Exponential function to use. Default to UFL exponential. :type exp: Callable :returns: Approximated Fermi-Dirac integral of order 1/2 at ``x``. :rtype: Any .. py:function:: d_density(x1: Any, x2: Any, Nmax: float, exp: collections.abc.Callable = ufl.exp) -> Any Neutral or ionized defect density function. :param x1: Reduced chemical potential of neutral or ionized defects. :type x1: Any :param x2: Reduced chemical potential of ionized or neutral defects. :type x2: Any :param Nmax: Maximum defect density. :type Nmax: float :param exp: Exponential function to use. Default to UFL exponential. :type exp: Callable :returns: Defect density. :rtype: Any .. note:: Pass respectively reduced chemical potentials of neutral and ionized defects to ``x1`` and ``x2`` to get neutral defect density, and vice versa to get ionized defect density. .. py:function:: ufl_mat2voigt4strain(eps: ufl.core.expr.Expr) -> ufl.core.expr.Expr Convert a UFL strain matrix to Voigt notation (a UFL vector). :param eps: Input strain matrix. :type eps: Expr :returns: Strain in Voigt notation. :rtype: Expr .. py:function:: young_poisson2stiffness(E: float, nu: float, dim: int) -> list[list[float]] Convert Young's modulus and Poisson's ratio to stiffness tensor in Voigt notation. :param E: Young's modulus. :type E: float :param nu: Poisson's ratio. :type nu: float :param dim: Dimension (1, 2, or 3). :type dim: int :returns: Stiffness tensor in Voigt notation. :rtype: list[list[float]] .. py:function:: ufl_tr_voigt(a: ufl.core.expr.Expr) -> ufl.core.expr.Expr Calculate the trace of a UFL matrix in Voigt notation. :param a: Input matrix in Voigt notation (a vector). :type a: Expr :returns: Trace of the input matrix. :rtype: Expr .. py:class:: L2RError(ref: dict[str, dolfinx.fem.Function], u: dict[str, dolfinx.fem.Function], jit_options: Optional[dict] = None, form_compiler_options: Optional[dict] = None, metadata: Optional[dict] = None) Class for computing L2 relative error between two sets of DOLFINx Functions. .. py:method:: setup(ref: dict[str, dolfinx.fem.Function], u: dict[str, dolfinx.fem.Function], jit_options: Optional[dict] = None, form_compiler_options: Optional[dict] = None, metadata: Optional[dict] = None) Compile forms for computing L2 relative error. :param ref: Reference DOLFINx Functions. :type ref: dict[str, fem.Function] :param u: DOLFINx Functions to compare against the reference. Only takes into account the items with the same keys as in ``ref``. :type u: dict[str, fem.Function] :param jit_options: Options to pass to just in time compiler :type jit_options: Optional[dict] :param form_compiler_options: Options to pass to the form compiler :type form_compiler_options: Optional[dict] :param metadata: Data to pass to the integration measure :type metadata: Optional[dict] .. py:method:: compute(eps: float = 1e-10) -> float Compute the L2 relative error between the two attached sets of DOLFINx Functions. :param eps: Small value to avoid division by zero. :type eps: float :returns: L2 relative error. :rtype: float .. py:class:: Projector(space: dolfinx.fem.FunctionSpace, petsc_options: Optional[dict] = None, jit_options: Optional[dict] = None, form_compiler_options: Optional[dict] = None, metadata: Optional[dict] = None) Projector for a given function. Solves Ax=b, where .. highlight:: python .. code-block:: python u, v = ufl.TrialFunction(Space), ufl.TestFunction(space) dx = ufl.Measure("dx", metadata=metadata) A = inner(u, v) * dx b = inner(function, v) * dx(metadata=metadata) :param function: UFL expression of function to project :param space: Space to project function into :param petsc_options: Options to pass to PETSc :param jit_options: Options to pass to just in time compiler :param form_compiler_options: Options to pass to the form compiler :param metadata: Data to pass to the integration measure .. attention:: This class is copied from http://jsdokken.com/FEniCS23-tutorial/src/approximations.html .. py:method:: reassemble_lhs() .. py:method:: assemble_rhs(h: ufl.core.expr.Expr) Assemble the right hand side of the problem .. py:method:: project(h: ufl.core.expr.Expr) -> dolfinx.fem.Function Compute projection using a PETSc KSP solver .. py:method:: __del__() .. py:class:: TMat Assemble and store temperature matrix for a given function space. It provides Cholesky decomposition and lumped temperature/mass matrix for calculating properly correlated noise out of an uncorrelated noise. The temperature matrix T_ij is defined as: ∫ T(x) v_i(x) v_j(x) dx, where T(x) is the temperature field and v_i(x), v_j(x) are the basis functions of the function space. :var lumpedT: The lumped temperature matrix (a vector). :vartype lumpedT: PETSc.Vec :var lumpedM: The lumped mass matrix (a vector). :vartype lumpedM: PETSc.Vec .. py:attribute:: lumpedT :type: petsc4py.PETSc.Vec .. py:attribute:: lumpedM :type: petsc4py.PETSc.Vec .. py:method:: setT(T: dolfinx.fem.Function | tuple[ufl.core.terminal.FormArgument | float, dolfinx.fem.FunctionSpace], jit_options: Optional[dict] = None, form_compiler_options: Optional[dict] = None, metadata: Optional[dict] = None) Attach the temperature field and assemble all the internal data for the first time. :param T: Temperature field, represented either as a :class:`dolfinx.fem.Function` or a tuple ``(T_, V)``, where ``T_`` is a UFL expression or a float (for constant temperature) and ``V`` is its (collapsed) function space. :type T: fem.Function | tuple[FormArgument | float, fem.FunctionSpace] :param jit_options: Options to pass to just-in-time compiler :type jit_options: Optional[dict] :param form_compiler_options: Options to pass to the form compiler :type form_compiler_options: Optional[dict] :param metadata: Data to pass to the integration measure :type metadata: Optional[dict] .. py:method:: factorT(type: str) Assemble and factor the temperature matrix for the attached temperature field with the current value. :param type: Type of factorization. E.g., "cholesky". :type type: str .. py:method:: solve_backwardT(b: petsc4py.PETSc.Vec, x: petsc4py.PETSc.Vec) Given the factored temperature matrix T = L U, solve the system U x = b. :param b: Right-hand side vector. :type b: PETSc.Vec :param x: Solution vector. :type x: PETSc.Vec .. caution:: Currently, parallel backward solving for Cholesky factorization through PETSc is only supported by MKL CPARDISO solver. .. py:method:: mat_solveM(B: petsc4py.PETSc.Mat, X: petsc4py.PETSc.Mat) Solve the system M X = B, where M is the mass matrix. :param B: Right-hand side matrix. :type B: PETSc.Mat :param X: Solution matrix. :type X: PETSc.Mat .. py:method:: assemble_lumpedT() Assemble the lumped temperature matrix :py:attr:`lumpedT` for the attached temperature field with the current value. .. hint:: No need to call this function after :py:meth:`~medypt.utils.TMat.setT` unless the temperature field has changed. .. py:method:: __del__() .. py:function:: gen_therm_noise(Tmat: TMat, dissipUmat: Any, rng: numpy.random.Generator, noise: dolfinx.fem.Function) Generate thermal noise based on the fluctuation-dissipation theorem. It is calculated as: M w = √2 A η L^T, where M is the mass matrix and A, L are the Cholesky factors of the temperature and dissipation matrices, respectively. η is a matrix of shape ``(Tmat_dim, dissipUmat_dim)`` containing uncorrelated standard normal random numbers, and w of the same shape is the *DOFs* of the proper thermal noise function. .. attention:: Currently only use method of lumped temperature and mass matrices for efficiency. :param Tmat: Assembled temperature matrix. :type Tmat: TMat :param dissipUmat: Upper-triangular Cholesky factor of the dissipation matrix. Can be a scalar for isotropic dissipation. :type dissipUmat: Any :param rng: Numpy random number generator. :type rng: Generator :param noise: Output noise field. It is the genuine thermal noise multiplied by √dt so that it is independent of dt, where dt is the time step size. :type noise: fem.Function