from fenics import * mesh = UnitSquareMesh(4, 4)
V = VectorFunctionSpace(mesh, "CG", 1, dim=2) P = FunctionSpace(mesh, "CG", 1) Up = VectorFunctionSpace(mesh, "R", 0, dim=2) Lp = VectorFunctionSpace(mesh, "R", 0, dim=2) S = MixedFunctionSpace([V, P, Up, Lp]) s0 = interpolate(Constant((0.0, 0.0, 0.0, 1.0, 2.0, 0.0, 0.0)), S) # extract a copy vp = Function(Up) assign(vp, s0.sub(2)) # doing some check # ... # assign a new value vp.vector()[0] = 4.0 vp.vector()[1] = 7.0 assign(s0.sub(2), vp) Regards, Simone On 26 Sep 2014, at 22:09, Debanjan Mukherjee <[email protected]> wrote: > Hello, > > I am working with a mixed function space created as follows: > > V = VectorFunctionSpace(mesh, "CG", 1, dim=2) > P = FunctionSpace(mesh, "CG", 1) > Up = VectorFunctionSpace(mesh, "R", 0, dim=2) > Lp = VectorFunctionSpace(mesh, "R", 0, dim=2) > > S = MixedFunctionSpace([V, P, Up, Lp]) > > Now, at one step in my calculation, I want to extract the vector from Up and > assign this into a numpy array to perform some additional analysis (which > involves a conditional check), and then set the values in a modified numpy > array and assign it to the vector in Up again. My attempt at doing this was > as follows: > > s0 = Function(S) > (u,p,vp,lp) = s.split() > vpNpy = vp.vector().array() > > and then modify the numpy array vpNpy into vpNpyNew and assign it back again > as follows: > > s1 = Function(S) > (u1,p1,vp1,lp1) = assign(s1.sub(2), vpNpyNew) > > This does not work, and I get an error message saying: > > "Cannot access a non-const vector from a subfunction" > > I tried looking through the Questions Forum - but I did not find a good > solution to this. > > Any help in this regard will be greatly appreciated. > > Regards > -- > Dr. Debanjan Mukherjee > > _______________________________________________ > fenics mailing list > [email protected] > http://fenicsproject.org/mailman/listinfo/fenics
_______________________________________________ fenics mailing list [email protected] http://fenicsproject.org/mailman/listinfo/fenics
