Subject: Re: [FFC-dev] [DOLFIN-dev] [Fwd: Re: [FEniCS-users] [Fwd: Re: tensor function]]

From: Kristian Oelgaard <[email-masked]>
Date: Thu, 01 Oct 2009 09:45:10 +0200

Quoting Patrick Riesen <[email-masked]>:

> i suppose the problem lies not in VTKFile.cpp.
>
> you consider a tensor function to have a rank >= 2,
> but for my tensor function the rank returned by e.g.
>
> e.function_space().element().value_rank()
>
> is 1. So it is rather a vector with 4 components having rank 1.
> possibly value_rank() in the FiniteElement class should return 2 even
> for a 2d tensor function, then your code in VTKFile is fine, but i
> didn't look into the FiniteElement classes.

The problem is in ffc/fem/createelement.py where a tensor element is treated as
a MixedElement which has rank 1. A simple workaround might be to attach
len(ufl_element._shape) as the rank of a mixed element and use this value when
we generate code for the value_rank() function.
What do you think Anders?

Kristian

> patrick
>
>
> -------- Original Message --------
> Subject: Re: [FEniCS-users] [Fwd: Re:  tensor function]
> Date: Wed, 30 Sep 2009 11:30:45 +0200
> From: Patrick Riesen <[email-masked]>
> CC: fenics <[email-masked]>
> References: <[email-masked]>
> <[email-masked]> <[email-masked]>
>
> seems that the tensor function for a 2d problem doesn't get to the
> output structure for a tensor function but to that for a vector function
> because it only has rank = 1. the write out for tensor is only invoked
> for a (3x3) tensor with rank = 2.
>
> Garth N. Wells wrote:
> >
> >
> > Patrick Riesen wrote:
> >> Hello,
> >>
> >> I don't know if this is a problem (dolfin0.9.2):
> >> When i was writing the TensorFunction e to a vtk-file after a solve
> >> by just using
> >>
> >> file << e
> >>
> >> the results i got were absolutely nonsense. I got three components,
> >> but none of them made any sense. This is why i was thinking i'm doing
> >> something completely wrong in my equations/code,
> >>
> >
> > VTK output of tensor functions is untested. I would be very surprised if
> > it did work. Dig around in VTKFile.cpp to see if you can find the problem.
> >
> > Garth
>
>
>
>
>
> >
> >> But if i first separate all components and write them to separate
> >> files as
> >>      Function e11 = e[0];
> >>      Function e12 = e[1];
> >>      Function e21 = e[2];
> >>      Function e22 = e[3];
> >>
> >>      fileE11 << e11;
> >>      fileE12 << e12;
> >>      fileE21 << e21;
> >>      fileE22 << e22;
> >>
> >> the solutions look reasonable (uff!..) and e12 and e21 are the same.
> >>
> >> regards,
> >> patrick
> >>
> >>
> >> Patrick Riesen wrote:
> >>> sorry, something went wrong, this would be the message:
> >>>
> >>> -------- Original Message --------
> >>> Subject: Re: [FEniCS-users] tensor function
> >>> Date: Wed, 30 Sep 2009 10:52:02 +0200
> >>> From: Patrick Riesen <[email-masked]>
> >>> To: Anders Logg <[email-masked]>
> >>> References: <[email-masked]>
> >>> <20090928122033.GH2528@olorin>    <[email-masked]>
> >>> <20090928140134.GA24982@olorin> <[email-masked]>
> >>>
> >>> or maybe not..
> >>>
> >>> i can not do
> >>>
> >>> Function e11 = e[0][0];
> >>>
> >>> i get
> >>>
> >>> error: no match for 'operator[]' in
> >>> 'dolfin::Function::operator[](dolfin::uint)(0u)[0]'
> >>>
> >>>
> >>> but i can get components doing
> >>>
> >>> e11 = e[0];
> >>> e12 = e[1];
> >>> e21 = e[2];
> >>> e22 = e[3];
> >>>
> >>>
> >>> i suppose for dolfin 0.9.3 it is the case as described below and for
> >>> dolfin 0.9.2 the case as above therefore ?
> >>> (which is what i am still using currently)
> >>>
> >>> greetings,
> >>> patrick
> >>>
> >>> Patrick Riesen wrote:
> >>>> this helped to sort things out, thanks again.
> >>>> patrick
> >>>>
> >>>> Anders Logg wrote:
> >>>>> On Mon, Sep 28, 2009 at 02:34:30PM +0200, Patrick Riesen wrote:
> >>>>>> Hi,
> >>>>>>
> >>>>>> thans for the quick response,
> >>>>>>
> >>>>>> therefore when V is the function space of the (e, v, p) element and
> >>>>>>
> >>>>>>> e[0] = (a, b)
> >>>>>>> e[1] = (c, d)
> >>>>>>> e[0][0] = a
> >>>>>>> e[0][1] = b
> >>>>>>> e[1][0] = c
> >>>>>>> e[1][1] = d
> >>>>>> therefore doing
> >>>>>>
> >>>>>> fs = new SubSpace(V,0,0)
> >>>>>>
> >>>>>> corresponds to the subspace for the two components (a, b) of the
> >>>>>> element  and not just for a, which then would be SubSpace(V,0,0,0) ?
> >>>>> Yes.
> >>>>>
> >>>>> --
> >>>>> Anders
> >>>> _______________________________________________
> >>>> FEniCS-users mailing list
> >>>> [email-masked]
> >>>> http://fenics.org/mailman/listinfo/fenics-users
> >>>
> >>> _______________________________________________
> >>> FEniCS-users mailing list
> >>> [email-masked]
> >>> http://fenics.org/mailman/listinfo/fenics-users
> >>
> >> _______________________________________________
> >> FEniCS-users mailing list
> >> [email-masked]
> >> http://fenics.org/mailman/listinfo/fenics-users
> >
> >
>
> _______________________________________________
> FEniCS-users mailing list
> [email-masked]
> http://fenics.org/mailman/listinfo/fenics-users
> _______________________________________________
> DOLFIN-dev mailing list
> [email-masked]
> http://www.fenics.org/mailman/listinfo/dolfin-dev
>

---

From: Anders Logg <[email-masked]>
Date: Thu, 1 Oct 2009 12:05:34 +0200

> What do you think Anders?

Feel free to try. Tensor functions are not very well supported at the
moment anyway. This is something we'll need to look into at some point
in the not so distant future. In particular, we're working on a solver
for the bidomain equations where we will need tensor fields.

It's not something I have time to look at just now though.

--
Anders

---

From: Kristian Oelgaard <[email-masked]>
Date: Thu, 01 Oct 2009 11:24:02 +0200

> It's not something I have time to look at just now though.

Neither do I, but I just added some code that should give the correct
value_rank() for tensor elements.
Patrick, could you try this out and report any problems related to vtk output?

Kristian

---

From: Patrick Riesen <[email-masked]>
Date: Thu, 1 Oct 2009 12:33:34 +0200

> Patrick, could you try this out and report any problems related to vtk output?

yes i will try it out tonight when i have my mac where i have the 
current development versions of fenics, i just have a fixed release of 
fenics09 on my institute workstation

thanks,
patrick

> FFC-dev mailing list
> [email-masked]
> http://www.fenics.org/mailman/listinfo/ffc-dev