medypt.model ============ .. py:module:: medypt.model .. autoapi-nested-parse:: Base class for phase-field models. Classes ------- .. autoapisummary:: medypt.model.ModelBase Module Contents --------------- .. py:class:: ModelBase Base class for phase-field models. .. py:attribute:: opts :type: dict[str, Any] A dictionary of numerical options. Initialized with default values and can be modified in subclasses. Contents include: * ``'has_thermal_noise'`` (bool): Whether to include thermal noise. Default is ``False``. * ``'rand_seed'`` (int): Random seed for thermal noise generation. Default is ``8347142``. * ``'quadr_deg'`` (int): Quadrature degree for numerical integration. Default is ``6``. * ``'petsc'`` (dict): A dictionary of PETSc solver options (see `this`_ for example). Default to use Newton nonlinear solver with MUMPS direct linear solver. * ``'t_step_rtol'`` (float): Relative tolerance for time discretization error. Default is ``0.01``. * ``'dt_min_rescalar'`` (float): Minimum factor to reduce time step upon failure. Default is ``0.2``. * ``'dt_max_rescalar'`` (float): Maximum factor to increase time step upon success. Default is ``4.0``. * ``'dt_reducer'`` (float): Factor to reduce time step rescalar. Default is ``0.9``. * ``'max_successive_fail'`` (int): Maximum number of successive failures before stopping. Default is ``100``. * ``'min_dt'`` (float): Minimum time step size. Default is ``1e-9``. * ``'max_dt'`` (float): Maximum time step size. Default is ``10.0``. * ``'save_period'`` (int): Time step period for saving solution. Default is ``1``. * ``'log_file_name'`` (str): File name for logging evolution. Default is ``evolution.txt``. * ``'sol_file_name'`` (str): File name for saving solution. Default is ``solution.xdmf``. If the suffix is not ``.xdmf`` or if no suffix, save solutions using :class:`dolfinx.io.VTXWriter` into a folder with the given name. * ``'verbose'`` (bool): Whether to print verbose output. Default is ``False``. .. _this: https://jsdokken.com/dolfinx-tutorial/chapter2/nonlinpoisson_code.html#newtons-method .. py:attribute:: params :type: dict[str, Any] A dictionary of physical parameters. Initialized to empty dictionary and should be set in subclasses. .. py:attribute:: mesh_data :type: dolfinx.io.gmsh.MeshData :class:`dolfinx.io.gmsh.MeshData` object containing mesh and boundary tags. Initialized to ``None`` and should be set in subclasses. .. py:attribute:: fields :type: dict[str, dolfinx.fem.Function] A dictionary mapping field names to their :class:`dolfinx.fem.Function` objects. Initialized to empty dictionary and should be set in subclasses. .. py:method:: load_mesh(comm: mpi4py.MPI.Comm, mesh: gmsh.model | str, **kwargs: Any) Load the mesh data. :param comm: MPI communicator. :type comm: MPI.Comm :param mesh: Gmsh model or a file name with the ``.msh`` or ``.xdmf`` format. :type mesh: gmsh.model | str :param kwargs: Additional keyword arguments passed to :py:func:`~medypt.utils.create_mesh`: * mesh_dim (int): Geometric dimension of the mesh. Required when ``mesh`` is a 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 when ``mesh`` is a gmsh model or a ``.msh`` file. * mesh_name (str): Name (identifier) of the mesh to read from the ``.xdmf`` file. Required when ``mesh`` is a ``.xdmf`` file. :type kwargs: dict[str, Any] :returns: ``None`` .. py:method:: set_bcs(bcs: list[tuple[str, int | collections.abc.Callable[[Any], Any], dolfinx.fem.Constant | dolfinx.fem.Function | numpy.ndarray | collections.abc.Callable[[Any], Any]]]) Set essential and natural boundary conditions. :param bcs: A list of tuples specifying boundary conditions. Each tuple is ``(name, b, bc)``, where: * ``name`` (str): The name of the field to which the boundary condition applies. * ``b`` (int | Callable[[Any], Any]): The tag of or a callable defining the boundary. The tag is only for boundaries with codimension 1. The callable is only for Dirichlet boundary conditions and will have a calling signature ``b(x)``, where ``x`` is the spatial coordinate treated as a 2D numpy array with shape (geom_dim, num_points). It should return a boolean array of shape ``(num_points,)`` indicating whether each point is on the boundary. * ``bc`` (Constant | Function | ndarray | Callable[[Any], Any]): A variable or callable defining the boundary condition. A variable is for Dirichlet boundary conditions, while a callable is for natural boundary conditions. The callable will have a calling signature ``bc(fields)``, where ``fields`` is treated as a tuple containing the splitted fields. It covers boundary conditions of the Neumann, Robin, and other general types implemented using the Nistche's trick [1]_ [2]_. The default is the zero-flux boundary condition. :type bcs: list[tuple[str, int | Callable[[Any], Any], fem.Constant | fem.Function | np.ndarray | Callable[[Any], Any]]] .. [1] https://doi.org/10.1016/j.camwa.2022.11.025 .. [2] https://doi.org/10.1103/PhysRevApplied.17.014042 .. tip:: One can make an evolving boundary condition by defining a global :py:class:`dolfinx.fem.Constant` or :py:class:`dolfinx.fem.Function` object and using it in the boundary-condition expression. One will then need to define an update function to update the value of the global object for a given time and pass that function to :py:meth:`~medypt.model.ModelBase.solve`. .. py:method:: create_problem() Create the finite element problem that is ready to be solved. .. py:method:: solve(tspan: float, dt: float | None = None, ics: dict[str, collections.abc.Callable[[Any], Any] | dolfinx.fem.Function] | None = None, update: collections.abc.Callable[[float], None] = lambda t: None) -> bool Solve the time-dependent finite-element problem. It uses adapative backward Euler time stepping algorithm for time integration. The time step is automatically controlled by the time discretization error specified in :py:attr:`~medypt.model.ModelBase.opts`. The advantage of backward Euler algorithm is its stability. :param tspan: Time span to solve over. :type tspan: float :param dt: Initial time step. Defaults to ``None`` to use the current time step. :type dt: float | None :param ics: Initial conditions, defined as a dictionary mapping field names to callables that return the initial value for that field. The callables have a calling signature ``ic(x)``, where ``x`` is the spatial coordinate treated as a 2D numpy array with shape ``(geom_dim, num_points)``. Defaults to ``None`` to use the current values of the fields. :type ics: dict[str, Callable[[Any], Any] | fem.Function] | None :param update: A callable to update any time-dependent parameters (must be defined a priori as global :py:class:`dolfinx.fem.Constant` or :py:class:`dolfinx.fem.Function`) at each time step. It takes the current time as the only argument. Defaults to a no-op function. :type update: Callable[[float], None] :returns: ``True`` if the solve completed successfully, ``False`` otherwise. :rtype: bool