Problem with test space for solving elasticity problem

dolfinx
1

Hi everyone,

I’m trying to solve some 2D linear elasticity inverse problem and have encountered with an issue I haven’t been able to work around yet. I have a formulation that uses vectorial test functions, but what I’m trying to compute is a 1D function. However, I’m not sure how to set this with a non-linear solver so I can compute this unknown function.

The minimal code goes as follow:

comm = MPI.COMM_WORLD

domain = dolfinx.mesh.create_rectangle(comm,
                            points=((0.0, 0.0), (30, 20)), n=(150, 100),
                            cell_type=mesh.CellType.triangle)

# create function spaces
vector_cg1 = ufl.VectorElement("CG", domain.ufl_cell(), 1)
V = fem.FunctionSpace(domain, vector_cg1) # for test space 

scalar_cg1 = ufl.FiniteElement("CG", domain.ufl_cell(), 1)
V_tilde = fem.FunctionSpace(domain, scalar_cg1)

v = ufl.TestFunction(V)
mu_fun = fem.Function(V_tilde)

# u_sol = some data function
# bc = some boundary condition for mu_fun

def epsilon(w):
  return (ufl.nabla_grad(w) + ufl.nabla_grad(w).T)/2

F = mu_fun*(ufl.inner(epsilon(u_sol),epsilon(v)) + ufl.nabla_div(u_sol)*ufl.nabla_div(v))*ufl.dx

problem = fem.petsc.NonlinearProblem(F, mu_fun, bcs=[bc])
solver = dolfinx.nls.petsc.NewtonSolver(comm, problem=problem)

solver.rtol = 1e-6
solver.convergence_criterion = "incremental"
solver.solve(mu_fun)

---------------------------------------------------------------------------
RuntimeError                              Traceback (most recent call last)
<ipython-input-55-1eea49d55584> in <module>
----> 1 num_its, converged = solver.solve(mu_fun)
      2 if converged:
      3     print(f" converged in {num_its} iterations")

2 frames
/usr/local/lib/python3.8/dist-packages/dolfinx/fem/assemble.py in set_bc(b, bcs, x0, scale)
    316 
    317     """
--> 318     _cpp.fem.set_bc(b, bcs, x0, scale)

RuntimeError: Size mismatch between b and x0 vectors.

I imagine this error comes from the dimension difference between the function spaces I’m using. Is this kind of formulation supported? Any hint or preexisting demo could be of great help.

Thanks in advance!

2

Please state the governing equation that you are trying to solve in Latex syntax. The variational form you have written here does not make sense.

3

Find \mu\in L^2(\Omega,\mathbb{R}) such that: \displaystyle\int_\Omega \mu(\varepsilon(u_h):\varepsilon(v)+(\nabla \cdot u_h)(\nabla \cdot v) )dx=0,\hspace{0.5cm} v\in H_0^1(\Omega,\mathbb{R}^2)
Where \varepsilon(u)=(\nabla u + (\nabla u)^T)/2 and \varepsilon(u):\varepsilon(v)=\displaystyle\sum_{ij}\varepsilon (u)_{ij}\varepsilon (v)_{ij}

4

This is not how you would solve an inverse problem. You need to formulate it as an optimization problem, ie.

\min_{\mu} \int_\Omega (u-u_{data})\cdot (u-u_{data})~\mathrm{d}x

such that

This can be solved by using the adjoint method.

5

Getting into more detail, to reduce time complexity for finding the problem solution, I was trying to turn it into some linear system. Given the basis from the finite spaces \{\phi_i\}_{i}^{dim V_h}\subseteq V_h\subseteq H_0^1(\Omega,\mathbb{R}^2), \{\psi_j\}_j^{dim U_h} \subseteq U_h\subseteq L^2(\Omega, \mathbb{R}), I’m trying to assemble the matrix: \mathcal{A}_{ij}=\displaystyle\int_\Omega\psi_j(\varepsilon(u_h):\varepsilon(\phi_i))+(\nabla\cdot u_h)(\nabla\cdot \phi_i)dx and solve: \mathcal{A}\mu=\vec{0} by minimizing ||\mathcal{A}^T\mathcal{A}\mu||^2_2. Can this be done via assembling the corresponding matrix and using PETSc?

6

Your problem doesn’t make sense for solution by Newton’s method. Perhaps just assemble your rectangular matrix \mathcal{A} and then do with that what you wish.

If you were to solve it as a PDE constrained minimisation problem as suggested by @dokken. That would assemble the \mathcal{A} \mathcal{A}^\top type problem you’re suggesting; and be well posed such that solution by Newton’s method should be straightforward.