Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Improve FEM API #457

Open
campospinto opened this issue Feb 4, 2025 · 0 comments
Open

Improve FEM API #457

campospinto opened this issue Feb 4, 2025 · 0 comments
Assignees
Labels
documentation enhancement New feature or request FEEC Finite element exterior calculus FEM API Objects describing finite element concepts multi-patch Next Release Must be in next release

Comments

@campospinto
Copy link
Collaborator

campospinto commented Feb 4, 2025

When trying to improve the compatibility between single and multi-patch simulations, I noticed a few issues with the Fem user interface. Below is a list, with some suggestions for improvement.

First a reminder of the different types of Fem spaces:

    - class psydac.fem.vector.ProductFemSpace:

      used for Vh = multi-patch FEM space
      
      then:
        
        - Vh.spaces[i] = FEM space on patch i
        - Vh.vector_space is of class psydac.linalg.block.BlockVectorSpace
          
    - class psydac.fem.vector.VectorFemSpace:
      
      used for Vh = single-patch, vector-valued FEM space
      
      then: 
      
        - Vh.spaces[i] = FEM space for (logical) component i
        - Vh.vector_space is of class psydac.linalg.block.BlockVectorSpace
        
    - class psydac.fem.tensor.TensorFemSpace:

      used for Vh = single-patch, scalar-valued tensor product FEM space 

      then:
      
        - Vh.spaces[i] = spline space along (logical) axis i
        - Vh.vector_space is of class psydac.linalg.stencil.StencilVectorSpace
      
    - class psydac.fem.splines.SplineSpace:

      used for Vh = scalar spline space on interval
      
      then: 
      
        - Vh has no attribute 'spaces' 
        - Vh.vector_space is of class 'psydac.linalg.stencil.StencilVectorSpace'

Then, a list of some problems that I saw:

  1. name ProductFemSpace is ambiguous: seems to be used only for multi-patch spaces, but a space of vector-valued fields
    (a VectorFemSpace) is also a product of FEM spaces...

    -> rename as MultipatchFemSpace ?

  2. "product" ambiguity is made worse by the flag function Vh.is_product which means either multi-patch or vector-valued
    (see for instance in the FemField constructor where its is understood as vector-valued...)

    -> replace by Vh.is_multipatch and Vh.is_vector_valued ?

  3. also the ProductFemSpace constructor may not return a ProductFemSpace, so the code is not easy to follow.

    -> use factory function instead ?

  4. Vh.spaces may return different types of spaces (as seen above) which makes the code not always easy to read, and hampers the single/multi-patch compatibility

    ->: add helper attributes? eg

    • Vh.patch_spaces[i]:
      if Vh is multi-patch space (ie ProductFemSpace), returns FEM space on patch i
      if Vh is single-patch (anything else) returns itself (and raise Error if i>0 ?)

    • Vh.component_spaces[i]:
      if Vh is multi-patch space (ie ProductFemSpace), raise error
      if Vh is single-patch, vector-valued (ie VectorFemSpace) returns space of component i
      if Vh is single-patch, scalar-valued (ie TensorFemSpace) returns itself (or error? only if i>0?)

    • Vh.axis_spaces[i]:
      if Vh is single-patch, scalar-valued tensor-product (ie TensorFemSpace), returns space for axis i
      else raise Error ?

  5. a similar problem occurs for a multipatch/vector-valued FemField uh (uh.fields may correspond to patch or component fields)

    -> add helper functions? eg

    • uh.patch_field(i):
      if uh is a multi-patch field returns FEM field on patch i
      else returns itself (raise error if i>0 ?)

    • Vh.component_field(i):
      needed ?

  6. several function or class docstrings indicate types of arguments which are wrong,
    eg:

    "L2 : SplineSpace" in class Projector_L2(GlobalProjector)

    "Vh : TensorFemSpace" in function construct_hcurl_conforming_projection
    (in feec.multipatch.non_matching_operators)

    -> better review of MR, use type hints, asserts ?

  7. the class DiscreteDerhamMultipatch is a subclass of DiscreteDerham (single patch by default)
    but it does not implement all the basic methods.
    In particular, the derivatives method which returns the differential operators is not implemented.

-> all the make sure that all the functions of DiscreteDerham are well implemented for DiscreteDerhamMultipatch.
In particular, for a multipatch dR sequence the derivatives should return the broken_derivatives (as operators, matrices).
should be understood as broken ones

  1. the DiffOperator class does not subclass FemLinearOperator, although they are very similar

-> Doing so will probably simplify a few things. Also,to_sparse_matrix (to convert a FemLinearOperator to sparse scipy format)
can probably be merged with tosparse

  1. the folder structure for single and multi-patch feec is not consistent:
    • DiscreteDerhamMultipatch in feec/multipatch/api
    • DiscreteDerham in api/feec

-> move the DiscreteDerhamMultipatch api in api/broken_feec ? Or api/feec/broken ?
(original idea was that classes defined in the api directory are exposed to the user)

There is also a related problem in Sympde:

  1. domain.interior.args returns objects of different types:

    • for a multi-patch domain: the subdomains as a tuple of logical domains
    • for a single patch domain: the name of the unique subdomain (the domain itself), as a tuple of strings,

    -> this is probably because domain.interior is a Union of subdomains.
    To fix this without changing how Unions are handled in Sympde, we could add a function/attribute domain.subdomains which always returns the subdomains as a tuple of logical domains

@campospinto campospinto changed the title FEM user interface not always clear improving the FEM api Feb 5, 2025
@campospinto campospinto self-assigned this Feb 5, 2025
@campospinto campospinto added enhancement New feature or request FEEC Finite element exterior calculus FEM API Objects describing finite element concepts documentation multi-patch labels Feb 5, 2025
yguclu added a commit that referenced this issue Mar 5, 2025
Update Psydac after the breaking change introduced by version 0.19.1 of
SymPDE: domain coordinates are real symbols.

Further, the new version is needed in order to solve #457.

---------

Co-authored-by: Elena Moral Sánchez <[email protected]>
@yguclu yguclu added the Next Release Must be in next release label Mar 11, 2025
@yguclu yguclu changed the title improving the FEM api Improve FEM API Mar 11, 2025
max-models added a commit to max-models/psydac-for-struphy that referenced this issue Mar 25, 2025
* fix bug gmres (pyccel#452)

Solves issue pyccel#451.

* CI Updates: Code coverage report, Ubuntu 24.04, PETSc 3.22 (pyccel#443)

Summary of changes:

- Integrate Codacy reporting
- Lock runner versions to `ubuntu-24.04` and `macos-14`. The
[`*-latest` runner
versions](https://github.com/actions/runner-images#available-images) are
moving targets; we wish to avoid surprise failing tests caused by runner
upgrades (e.g. see [macOS
issue](pyccel#421) from last year) .
- Upgrade PETSc to [v3.22](https://petsc.org/release/changes/322/).

* Expose inter-/histopolation matrices in GlobalProjector, add unit tests (pyccel#347)

Fixes pyccel#346.

The inter-/histopolation matrices are now exposed in
`GlobalProjector.imat_kronecker` as a `LinearOperator` object of type
`KroneckerStencilMatrix` (in the case where both domain and codomain are
scalar spaces) or `BlockLinearOperator` (in the case where at least one
space is vector valued).

The unit tests in `feec/tests/test_commuting_projections.py` have been
enhanced: they now also test how close `imat_kronecker` is to being the
right-inverse of the `solver`attribute (which uses a Kronecker product
of 1D linear solvers).

---------

Co-authored-by: Yaman Güçlü <[email protected]>
Co-authored-by: Max Lindqvist <[email protected]>
Co-authored-by: Martin Campos Pinto <[email protected]>

* Avoid installing h5py from cache (pyccel#430)

Add the flag `--no-cache-dir` to `pip` in order to guarantee
that `h5py` is always built during installation. See:

  - https://stackoverflow.com/questions/45594707/what-is-pips-no-cache-dir-good-for;
  - https://docs.h5py.org/en/stable/build.html#custom-installation.

Re-using a prebuilt `h5py` (*i.e.* cached `h5py`) could become an issue
when `h5py` dependencies upgrade. For example, If people have a cached
version of the h5py library which was compiled using NumPy < 2.0 and try
to install Psydac in a clean Python environment, they end up with a broken
version. This is because the cached version of h5py tries to call C functions
from NumPy 2.0 which are not ABI compatible.

---------

Co-authored-by: Yaman Güçlü <[email protected]>

* Allow keyword parameters in the dot method of `MatrixFreeLinearOperator` (pyccel#456)

Allow the creation of a `MatrixFreeLinearOperator` object with a dot
function that takes multiple keyword arguments, as for instance function
parameters.

* Small fixes for operators with complex (pyccel#455)

- Make `ScaledLinearOperator` work with a complex scalar;
- Take real part of `r.dot(r)` to avoid complex warning in `sqrt`;
- Add unit test.

* Fix Greville point collocation matrix for single grid point case (pyccel#446)

Previously, in the case when `self.greville.size == 1`,
`collocation_matrix` would return an empty matrix.

Now, if the grid has a single Greville point, the collocation matrix is
simplified to a 1x1 matrix with a value of 1.

* Use SymPDE version 0.19.1 (pyccel#445)

Update Psydac after the breaking change introduced by version 0.19.1 of
SymPDE: domain coordinates are real symbols.

Further, the new version is needed in order to solve pyccel#457.

---------

Co-authored-by: Elena Moral Sánchez <[email protected]>

* Use CI from psydac-for-struphy

* Install with .[test]

---------

Co-authored-by: Elena Moral Sánchez <[email protected]>
Co-authored-by: Paul Rigor <[email protected]>
Co-authored-by: Stefan Possanner <[email protected]>
Co-authored-by: Yaman Güçlü <[email protected]>
Co-authored-by: Martin Campos Pinto <[email protected]>
@yguclu yguclu pinned this issue Apr 2, 2025
spossann added a commit to max-models/psydac-for-struphy that referenced this issue Apr 4, 2025
* fix bug gmres (pyccel#452)

Solves issue pyccel#451.

* CI Updates: Code coverage report, Ubuntu 24.04, PETSc 3.22 (pyccel#443)

Summary of changes:

- Integrate Codacy reporting
- Lock runner versions to `ubuntu-24.04` and `macos-14`. The
[`*-latest` runner
versions](https://github.com/actions/runner-images#available-images) are
moving targets; we wish to avoid surprise failing tests caused by runner
upgrades (e.g. see [macOS
issue](pyccel#421) from last year) .
- Upgrade PETSc to [v3.22](https://petsc.org/release/changes/322/).

* Expose inter-/histopolation matrices in GlobalProjector, add unit tests (pyccel#347)

Fixes pyccel#346.

The inter-/histopolation matrices are now exposed in
`GlobalProjector.imat_kronecker` as a `LinearOperator` object of type
`KroneckerStencilMatrix` (in the case where both domain and codomain are
scalar spaces) or `BlockLinearOperator` (in the case where at least one
space is vector valued).

The unit tests in `feec/tests/test_commuting_projections.py` have been
enhanced: they now also test how close `imat_kronecker` is to being the
right-inverse of the `solver`attribute (which uses a Kronecker product
of 1D linear solvers).

---------

Co-authored-by: Yaman Güçlü <[email protected]>
Co-authored-by: Max Lindqvist <[email protected]>
Co-authored-by: Martin Campos Pinto <[email protected]>

* Avoid installing h5py from cache (pyccel#430)

Add the flag `--no-cache-dir` to `pip` in order to guarantee
that `h5py` is always built during installation. See:

  - https://stackoverflow.com/questions/45594707/what-is-pips-no-cache-dir-good-for;
  - https://docs.h5py.org/en/stable/build.html#custom-installation.

Re-using a prebuilt `h5py` (*i.e.* cached `h5py`) could become an issue
when `h5py` dependencies upgrade. For example, If people have a cached
version of the h5py library which was compiled using NumPy < 2.0 and try
to install Psydac in a clean Python environment, they end up with a broken
version. This is because the cached version of h5py tries to call C functions
from NumPy 2.0 which are not ABI compatible.

---------

Co-authored-by: Yaman Güçlü <[email protected]>

* Allow keyword parameters in the dot method of `MatrixFreeLinearOperator` (pyccel#456)

Allow the creation of a `MatrixFreeLinearOperator` object with a dot
function that takes multiple keyword arguments, as for instance function
parameters.

* Small fixes for operators with complex (pyccel#455)

- Make `ScaledLinearOperator` work with a complex scalar;
- Take real part of `r.dot(r)` to avoid complex warning in `sqrt`;
- Add unit test.

* Fix Greville point collocation matrix for single grid point case (pyccel#446)

Previously, in the case when `self.greville.size == 1`,
`collocation_matrix` would return an empty matrix.

Now, if the grid has a single Greville point, the collocation matrix is
simplified to a 1x1 matrix with a value of 1.

* Use SymPDE version 0.19.1 (pyccel#445)

Update Psydac after the breaking change introduced by version 0.19.1 of
SymPDE: domain coordinates are real symbols.

Further, the new version is needed in order to solve pyccel#457.

---------

Co-authored-by: Elena Moral Sánchez <[email protected]>

* Improve FEM API (names, properties, constructors...) & plot_field() (pyccel#468)

Improve the FEM API to make it more user-friendly and compatible between
single-patch and multi-patch APIs (fix pyccel#459 and fix pyccel#460):

* Rename abstract property `vector_space` of `FemSpace` as `coeff_space`
* Rename `ProductFemSpace` as `MultipatchFemSpace`
* Do not use class constructors as factories (previously
`ProductFemSpace()` could return a `VectorFemSpace` object)
* Provide more specific helper properties for FEM spaces, such as
- `is_vector_valued` and `is_multipatch` (also remove a method
`is_scalar` which was implemented in _some_ FEM spaces but not used)
- `patch_spaces`, `component_spaces` and `axis_spaces` which always
return tuples of Fem spaces (or an error)
* Provide more specific helper properties for FEM fields, such as
- `component_fields` and `patch_fields` and `axis_fields` which always
return tuples of FEM fields (or an error)
* Improve the `VectorFemSpace` class (in particular its constructor) to
better match the `FemSpace` abstract interface.

Since these changes are also motivated by a better compatibility between
single-patch and multi-patch APIs (see issue pyccel#331), this PR also
modifies the `plot_field` function (previously in
`psydac.feec.multipatch.plotting_utilities`), which is:

* made compatible with single patch and multi patch fields
* renamed as `plot_field_2d` since it only handles 2d fields
* tested (not the content but the runs) in single and multipatch
configurations for scalar and vector-valued FEM fields
* moved with its module `psydac.feec.multipatch.plotting_utilities` to
`psydac.fem.plotting_utilities` since it is not specific to FEEC or
multipatch.

---------

Co-authored-by: Yaman Güçlü <[email protected]>

* Enable mpi4py 4.0 (pyccel#478)

Closes pyccel#420.

Fix usage of function `Translate_ranks` which had
non-backward-compatible changes in mpi4py version 4.0:
https://mpi4py.readthedocs.io/en/stable/changes.html.

---------

Co-authored-by: Stefan Possanner <[email protected]>
Co-authored-by: Max Lindqvist <[email protected]>
Co-authored-by: Yaman Güçlü <[email protected]>

* Removed dependency <4 for mpi4py

* Set ubuntu-latest and macos-14 in pyproject.toml

* Removed parallel hdf5

* Install struphy from rename-vector_space-to_coeff_space in CI

* Removed python 3.13 from the CI

* 3.13 --> 3.12 for macos CI

* Removed Run MPI tests with Pytest

* Removed codacy

* Support Python 3.13 (pyccel#475)

Closes pyccel#476:

* Do not restrict the maximum Python version to 3.12 in `pyproject.toml`
* Require `Cython >= 3` to avoid `h5py` installation crash w/ Python 3.13
* Require `sympde == 0.19.2` which supports Python 3.13
* Run unit tests with Python 3.13 too

* Added Set up environment variables

* Commented out macos tests

* Install struphy from devel branch

---------

Co-authored-by: Elena Moral Sánchez <[email protected]>
Co-authored-by: Paul Rigor <[email protected]>
Co-authored-by: Stefan Possanner <[email protected]>
Co-authored-by: Yaman Güçlü <[email protected]>
Co-authored-by: Martin Campos Pinto <[email protected]>
Co-authored-by: Stefan Possanner <[email protected]>
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
documentation enhancement New feature or request FEEC Finite element exterior calculus FEM API Objects describing finite element concepts multi-patch Next Release Must be in next release
Projects
None yet
Development

No branches or pull requests

2 participants