This is a read only copy of the old FEniCS QA forum. Please visit the new QA forum to ask questions

store several solutions in one variable

0 votes

I want to calculate the Gramian matrix of several solution, in order to do that I need to store my solution in one variable so after calculating all the solutions i can load/call the variable. so it should be something like :
u = Function(V)
solve(a == L, u, bc)
B[i]=u
how can I define B ?

I know there is Time series class but I rather not to store u.vector() instead of u.

from __future__ import division
from dolfin import *
import pickle
import numpy, scipy.io
import csv
E0=[69000000000,65000000000];
Magnitude=-100;nu = 0.3
mesh=UnitSquareMesh(10, 10)
V = VectorFunctionSpace(mesh, "Lagrange", 1)
WW = VectorFunctionSpace(mesh, "Lagrange", 1, 2)
B = Function(WW)
for i in range(0,len(E0)):
     class Bottom(SubDomain):
              def inside(self, x, on_boundary):
                    return (near(x[1], 0.0) )
      class Top(SubDomain):
               def inside(self, x, on_boundary):
                     return (near(x[1], 1.0) )
     top = Top()
     bottom = Bottom()
     boundaries = FacetFunction("uint", mesh)
     boundaries.set_all(0)
     top.mark(boundaries, 2)
     bottom.mark(boundaries, 4)
     bc = DirichletBC(V, (0.0, 0.0), boundaries,4)
     ds = Measure("ds")[boundaries]
     u = TrialFunction(V)
     v = TestFunction(V)
     f = Expression(("0.0", "scale"),  scale = Magnitude )
     E=E0[i]
    mu = E / (2.0*(1.0 + nu))
    lmbda = E*nu / ((1.0 + nu)*(1.0 - 2.0*nu))
   def sigma(u):
         return 2.0*mu*sym(grad(u)) + lmbda*tr(sym(grad(u)))*Identity(v.cell().d)
   a = inner(sigma(u), sym(grad(v)))*dx
   L = inner(f,v)*ds(2)
   u = Function(V)
   solve(a == L, u, bc)
   assign(B.sub(i), u)

so for just sake of my second question if I just store each solution in different variable and then want to calculate the covariance matrix of my solutions, if i want to use the same approach as "http://fenicsproject.org/qa/1826/compute-gramian-matrix" which means using assemble(U1U2dx) I get error. and when i am using a12=inner(U1,U2)*dx, a12 is ufl form but I expect that have scalar value for a12

from __future__ import division
from dolfin import *
import pickle
import numpy, scipy.io
import csv
E0=[69000000000,65000000000];
Magnitude=-100;nu = 0.3
mesh=UnitSquareMesh(10, 10)
V = VectorFunctionSpace(mesh, "Lagrange", 1)
WW = VectorFunctionSpace(mesh, "Lagrange", 1, 2)
B = Function(WW)
for i in range(0,len(E0)):
     class Bottom(SubDomain):
               def inside(self, x, on_boundary):
                      return (near(x[1], 0.0) )

      class Top(SubDomain):
              def inside(self, x, on_boundary):
                    return (near(x[1], 1.0) )
      top = Top()
      bottom = Bottom()
      boundaries = FacetFunction("uint", mesh)
      boundaries.set_all(0)
      top.mark(boundaries, 2)
      bottom.mark(boundaries, 4)
      bc = DirichletBC(V, (0.0, 0.0), boundaries,4)
      ds = Measure("ds")[boundaries]
       u = TrialFunction(V)
       v = TestFunction(V)
       f = Expression(("0.0", "scale"),  scale = Magnitude )
      E=E0[i]
      mu = E / (2.0*(1.0 + nu))
      lmbda = E*nu / ((1.0 + nu)*(1.0 - 2.0*nu))
       def sigma(u):
             return 2.0*mu*sym(grad(u)) + lmbda*tr(sym(grad(u)))*Identity(v.cell().d)
       a = inner(sigma(u), sym(grad(v)))*dx
       L = inner(f,v)*ds(2)
       u = Function(V)
       solve(a == L, u, bc)

       if i==0:
           U1=u
       else:
           U2=u

alpha = inner(U1,U2)*dx

Any hints, or ideas ?

asked Dec 18, 2013 by Bahram FEniCS Novice (400 points)
edited Dec 19, 2013 by Bahram

1 Answer

+1 vote

If the solutions are in the same space you can store them in a VectorFunction:
(requires development version of DOLFIN)

W = VectorFunctionSpace(mesh, space, deg, num_solutions)
B = Function(W)
for i, (a, L, bc) in enumerate(linear_problems):
    solve(a==L, u, bc)
    assign(B.sub(i), u)

You need to provide a minimal but runnable code for us to help you with the second question.

answered Dec 19, 2013 by johanhake FEniCS Expert (22,480 points)

Thank you very much for your help, I tried to store the solution in a way that you suggested but i got an error and I put the code in the main question. also for the second question I put the minimal running code in the main question.
Thank you again for your help.

...