← Back to team overview

ufl team mailing list archive

Re: [Question #215432]: True User-defined Operators for UFL

 

Question #215432 on UFL changed:
https://answers.launchpad.net/ufl/+question/215432

    Status: Answered => Open

David M. Rogers is still having a problem:
  Thanks for the quick reply.  I was thinking along those lines (and
would of course need dp/dc), but was unsure where to start.  Are you
saying that I can do the following to the CahnHilliard example?

c = Function(V) # to be solved for, index J
p = Coefficient(V) # represents q = p(c), indices I,J
v = TestFunction(V), index I

L = dot(p*nabla_grad(c), nabla_grad(v))*dx

ct = TrialFunction(V) # index K
dpdc = Coefficient(V) # represents s = dp/dc, indices I,J,K

# a = derivative(L, c, ct) # Gateaux derivative in dir. of ct
# From Sec. 1.2.4 of the Fenics Manual:
a = dot(p*nabla_grad(ct), nabla_grad(v))*dx \
  + dot(dpdc*ct*nabla_grad(x), nabla_grad(v))*dx

class UserFunction(Expression): ...
class UserDeriv(Expression): ...

fp = UserFunction(c) # stores location of c's data
fdpdc = UserDeriv(c) # stores location of c's data

# Class for interfacing with the Newton solver
class CahnHilliardEquation(NonlinearProblem):
    def __init__(self, a, L, p, dpdc):
        NonlinearProblem.__init__(self)
        self.L = L
        self.a = a
        self.p = p
        self.dpdc = dpdc
        self.reset_sparsity = True
    def F(self, b, x):
        self.p.interpolate(f) # sees new 'c' when called?
        assemble(self.L, tensor=b)
    def J(self, A, x):
        self.p.interpolate(fp)
        self.dpdc.interpolate(fdpdc)
        assemble(self.a, tensor=A, reset_sparsity=self.reset_sparsity)
        self.reset_sparsity = False

  That would probably make the best solution as long as JIT compiling
will not be re-done.  Since  I'm defining p(c) and dpdc as 'Coefficient'
fields, this is a Dolfin issue, so I'll send further questions their
way.  If it works I'll try and send the example back to the Fenics
project.

Thanks!
~ David

-- 
You received this question notification because you are a member of UFL
Team, which is an answer contact for UFL.