dolfin team mailing list archive

[Johan Hake <johan.hake@xxxxxxxxx>]

On Monday 25 January 2010 10:01:18 Garth N. Wells wrote:
> Johan Hake wrote:
> > On Monday 25 January 2010 03:10:46 Garth N. Wells wrote:
> >> Johan Hake wrote:
> >>> On Sunday 24 January 2010 16:39:14 Garth N. Wells wrote:
> >>>> Johan Hake wrote:
> >>>>> On Sunday 24 January 2010 16:12:37 Garth N. Wells wrote:
> >>>>>> Johan Hake wrote:
> >>>>>>> On Sunday 24 January 2010 00:03:41 Garth N. Wells wrote:
> >>>>>>>> Johan Hake wrote:
> >>>>>>>>> On Saturday 23 January 2010 14:55:08 Garth N. Wells wrote:
> >>>>>>>>>> Johan Hake wrote:
> >>>>>>>>>>> On Saturday 23 January 2010 08:42:14 Garth N. Wells wrote:
> >>>>>>>>>>>> Is it correct that behind the scenes that
> >>>>>>>>>>>>
> >>>>>>>>>>>>   U0 = Function(V)
> >>>>>>>>>>>>   U  = Function(V)
> >>>>>>>>>>>>   U0.vector()[:] = U.vector()[:]
> >>>>>>>>>>>>
> >>>>>>>>>>>> involves a GenericVector::get(..) call and a
> >>>>>>>>>>>> GenericVector::set(..) call? If so, it isn't ideal since it
> >>>>>>>>>>>> introduces unnecessary new/delete operations and unnecessary
> >>>>>>>>>>>> copying of data.
> >>>>>>>>>>>
> >>>>>>>>>>> None of GenericVector::get(..) or GenericVector::set(..) are
> >>>>>>>>>>> invoked, see __getslice__ and __setslice__ in la_post.i.
> >>>>>>>>>>>
> >>>>>>>>>>>    U0.vector()[:]
> >>>>>>>>>>>
> >>>>>>>>>>> involves
> >>>>>>>>>>>
> >>>>>>>>>>>   GenericVector::operator =(..)
> >>>>>>>>>>>
> >>>>>>>>>>> and
> >>>>>>>>>>>
> >>>>>>>>>>>   U.vector()[:]
> >>>>>>>>>>>
> >>>>>>>>>>> involves
> >>>>>>>>>>>
> >>>>>>>>>>>   GenericVector::copy()
> >>>>>>>>>>>
> >>>>>>>>>>> However the latter is unnecessary as you instead can do:
> >>>>>>>>>>>
> >>>>>>>>>>>   U0.vector()[:] = U.vector()
> >>>>>>>>>>>
> >>>>>>>>>>> invoking the assignment operator of U0's vector with U's
> >>>>>>>>>>> vector.
> >>>>>>>>>>
> >>>>>>>>>> What happens if I do
> >>>>>>>>>>
> >>>>>>>>>>   x = U.vector()[:]
> >>>>>>>>>
> >>>>>>>>> It just triggers the copy method of GenericVector, which is the
> >>>>>>>>> same behavior as for other itterable Python types.
> >>>>>>>>>
> >>>>>>>>>> ? Is x a numpy array?
> >>>>>>>>>
> >>>>>>>>> No you need to call array() to accomplish that.
> >>>>>>>>
> >>>>>>>> OK. What I'm trying to do is
> >>>>>>>>
> >>>>>>>>     # Get vectors
> >>>>>>>>     u_vec  = u.vector()[:]
> >>>>>>>>     u0_vec = u0.vector()[:]
> >>>>>>>>     v0_vec = v0.vector()[:]
> >>>>>>>>     a0_vec = a0.vector()[:]
> >>>>>>>
> >>>>>>> You should not need to make a copy of the vectors here.
> >>>>>>
> >>>>>> How can I avoid it?
> >>>>>
> >>>>> a_vec and v_vec are new vectors. None of the four vectors below get
> >>>>> modified by the a_vec and v_vec expressions so no need of copying,
> >>>>> and the v0 and a0 assignment should work with GenericVectors too.
> >>>>
> >>>> Do you mean that just
> >>>>
> >>>>   a_vec = 1.0/(2.0*beta)*((u - u0 - v0*dt)/(0.5*dt*dt) \
> >>>>            - (1.0-2.0*beta)*a0 )
> >>>>
> >>>> where u and u0 are GenericVectors should work?
> >>>
> >>> Have you tried?
> >>
> >> Can I somehow get a 'reference' to the vector so I don't have to use
> >> u.vector()[:] in the expressions?
> >
> > That is what u.vector() gives you.
> 
> I figured it out eventually. There's nothing quite like a '&' symbol . .
> . . .

Good!

In Python every name in a namespace is a reference(pointer) to an instance. 
Johan

> Garth
> 
> > Johan
> >
> >> Garth
> >>
> >>> As long as the rest (besides a0, which I assume also is a
> >>> GenericVector) are scalars everything should just work. The Python LA
> >>> interface (at least for GenericVector) should work more or less as the
> >>> NumPy interface which I think is nice :)
> >>>
> >>> We cannot take 1./v, where v is a GenericVector.
> >>>
> >>>>> Do you get any error messages?
> >>>>
> >>>> No errors. What I have now seems to work fine.
> >>>
> >>> Ok, and that is because what you do obviously works for NumPy arrays.
> >>>
> >>> Johan
> >>>
> >>>> Garth
>