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? 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 > > >> # Update acceleration and velocity > >> a_vec = (1.0/(2.0*beta))*( (u_vec - u0_vec - v0_vec*dt)/(0.5*dt*dt)\ > >> - (1.0-2.0*beta)*a0_vec ) > >> v_vec = dt*((1.0-gamma)*a0_vec + gamma*a_vec) + v0_vec > >> > >> > >> > >> Garth > >> > >>>> a_vec = 1.0/(2.0*beta)*((u_vec-u0_vec - v0_vec*dt)/(0.5*dt*dt) \ > >>>> - (1.0-2.0*beta)*a0_vec) > >>>> v_vec = dt*((1.0-gamma)*a0_vec + gamma*a_vec) + v0_vec > >>> > >>> The expression cause a lot of intermediate vectors, but should work, if > >>> beta, gamma and dt are all scalars. > >>> > >>>> # Update (U0 <- U0) > >>>> v0.vector()[:] = v_vec > >>>> a0.vector()[:] = a_vec > >>> > >>> This should work. > >>> > >>> Johan > >>> > >>>> I guess I need to call array() to make this work? > >>>> > >>>> Garth > >>>> > >>>>> Johan > >>>>> > >>>>>> Garth > >>>>>> > >>>>>>> Johan > _______________________________________________ Mailing list: https://launchpad.net/~dolfin Post to : [email protected] Unsubscribe : https://launchpad.net/~dolfin More help : https://help.launchpad.net/ListHelp
