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

Use of PWMatrixCoefficient to define matrix coefficients per attribute #253

Open
mdavids-cfs opened this issue Sep 16, 2024 · 3 comments
Open

Comments

@mdavids-cfs
Copy link

mdavids-cfs commented Sep 16, 2024

Hi all,
I have been trying to implement a PWMatrixCoefficient but I am having trouble with the interface (also couldn't find an example for this). The goal is to assign a MatrixConstantCoefficient per mesh attribute.
In normal MFEM (C++) this is done by creating an Array Array<MatrixCoefficient *> coefs(0); and then assigning the matrix coefficients. I thought that PyMFEM can interpret normal lists of tuples as list-equivalents, i.e., I tried something like this:

sigma_all_coefs = []
sigma_attr = mfem.intArray(max_attr)
for ti,tensor in enumerate(sigma_all): 
    # just for testing
    temp = mfem.DenseMatrix(dim)
    temp.Assign(0.0)

    # add matrix coefficient to list
    sigma_all_coefs.append(mfem.MatrixConstantCoefficient(temp)  )
    sigma_attr[ti] = ti+1   
    
# Create PW Matrix Coefficient
sigmaCoef = mfem.PWMatrixCoefficient(dim, sigma_attr, sigma_all_coefs, False)

which fails when calling the mfem.PWMatrixCoefficient constructor. Also tried another approach:

# Create PW Matrix Coefficient
sigmaCoef = mfem.PWMatrixCoefficient(dim, False)

for ti, tensor in enumerate(sigma_all):
    # Set coefficient
    sigmaCoef.UpdateCoefficient( ti+1, mfem.MatrixConstantCoefficient( mfem.DenseMatrix(tensor) ) )

which fails when attempting to assemble a bilinear operator, that utilizes the PWMatrixCoefficient.

Any thoughts on how this could be achieved?
Big thanks and warm regards,
Mathias

@v-dobrev
Copy link
Member

One option, as discussed in #57, may be to derive your own class from mfem.MatrixPyCoefficient as illustrated in https://github.com/mfem/PyMFEM/blob/master/examples/wrapper_example/matrix_coefficient.py.

Unfortunately, I don't know if one can create an equivalent of Array<MatrixCoefficient*> in python and I don't think sigma_all_coefs, as constructed in your first code snippet, has the right type. Also, in the second code snippet, when calling sigmaCoef.UpdateCoefficient the second argument in C++ is MatrixCoefficient & and I suspect that on the python side using mfem.MatrixConstantCoefficient is not correct but I don't know why -- I'm not familiar with how SWIG translates things.

@sshiraiwa
Copy link
Member

@mdavids-cfs, @v-dobrev This is going to be address in https://github.com/mfem/PyMFEM/tree/coefficient-arrray-dev, where MatrixCoefficientPtrArray is available. The second option may work in the current master, if you keep the argument
objects ( mfem.MatrixConstantCoefficient and mfem.DenseMatrix) from being gc-ed.

@mdavids-cfs
Copy link
Author

mdavids-cfs commented Oct 2, 2024

Hi @v-dobrev and @sshiraiwa
I should've responded earlier; I did in fact get it to working by storing MatrixConstantCoefficient and DenseMatrix in lists, so they don't get overwritten, and then using the second approach. This worked pretty well!
Also, big thanks for putting together this python interface to MFEM, it's really helpful!

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

3 participants