← Back to team overview

dolfin team mailing list archive

Re: [Bug 681415] Re: Problem using sub-spaces

 

On Tue, Nov 30, 2010 at 09:57:42AM -0000, Garth Wells wrote:
>
> On 29/11/10 22:51, Johan Hake wrote:
> > On Monday November 29 2010 14:39:33 Anders Logg wrote:
> >> On Mon, Nov 29, 2010 at 09:13:11PM -0000, Johan Hake wrote:
> >>> On Monday November 29 2010 12:59:55 Anders Logg wrote:
> >>>> On Mon, Nov 29, 2010 at 07:08:16PM -0000, Garth Wells wrote:
> >>>>> On 29/11/10 18:50, Johan Hake wrote:
> >>>>>> On Monday November 29 2010 10:23:43 Mikael Mortensen wrote:
> >>>>>>> On 29 November 2010 17:06, Johan Hake<681415@xxxxxxxxxxxxxxxxxx>
> >>>
> >>> wrote:
> >>>>>>>> On Monday November 29 2010 00:35:27 Mikael Mortensen wrote:
> >>>>>>>>> Please don't.
> >>>>>>>>>
> >>>>>>>>> The subspace functionality is used for example to compute the
> >>>>>>>>> streamfunction, where the subspace of the velocity
> >>>>>>>>> VectorFunctionSpace is used. I also use the subspace
> >>>>>>>>> functionality to assemble into parts of the Navier-Stokes
> >>>>>>>>> coefficient matrix componentwise. For example
> >>>>>>>>>
> >>>>>>>>> uc = TrialFunction(V.sub(0))
> >>>>>>>>> vc = TestFunction(V.sub(0))
> >>>>>>>>> a = assemble(inner(grad(uc), grad(vc))*dx)
> >>>>>>>>>
> >>>>>>>>> and then I can copy this 'small' matrix to the three diagonal
> >>>>>>>>> slots of the 3*3 Laplacian matrix that makes up for the whole
> >>>>>>>>> velocity vectorfield. I understand this is not the regular way
> >>>>>>>>> of assembling, but it is much faster than switching V.sub(0)
> >>>>>>>>> above to V and assembling the whole thing by itself.
> >>>>>>>>
> >>>>>>>> What hinders you to do:
> >>>>>>>>    Vs = FunctionSpace(mesh, "CG", 1)
> >>>>>>>>    V  = MixedFunctionSpace([Vs]*3)
> >>>>>>>>    ...
> >>>>>>>>    uc = TrialFunction(Vs)
> >>>>>>>>    vc = TestFunction(Vs)
> >>>>>>>>
> >>>>>>>>     a = assemble(inner(grad(uc), grad(vc))*dx)
> >>>>>>>
> >>>>>>> Nothing really. I used to do that. However, I want to be able to
> >>>>>>> grab the velocity subpace when I'm computing the streamfunction
> >>>>>>> as I want to use the same functionspace for the streamfunction as
> >>>>>>> a velocity-component to avoid any interpolation. I understand
> >>>>>>> that it's not strictly necessary, it's just a nice functionality
> >>>>>>> the way I see it and I can use the same
> >>>>>>> streamfunction-code for all velocity-functionspaces.
> >>>>>>
> >>>>>> Ok.
> >>>>>>
> >>>>>> What if we keep the original subspaces (without offset) in a
> >>>>>> protected list.
> >>>>>>
> >>>>>>     _spaces = [Vs0, Vs1, Vs2]
> >>>>>>
> >>>>>> Then we add a kwarg: 'with_offset=True' to the 'sub' method. If
> >>>>>> 'with_offset' is False we just return the "original" subspace
> >>>>>> without offset. You then need
> >>>>>>
> >>>>>> to change your code to:
> >>>>>>     uc = TrialFunction(V.sub(0, False))
> >>>>>>     vc = TestFunction(V.sub(0, False))
> >>>>>>     a = assemble(inner(grad(uc), grad(vc))*dx)
> >>>>>>
> >>>>>> Then we make it impossible to instantiate a Function,
> >>>>>> TrialFunction, TestFunction with a FunctionSpace including an
> >>>>>> offset.
> >>>>>>
> >>>>>> How do I read out from a FunctionSpace that it includes a ufc
> >>>>>> offset?
> >>>>>
> >>>>> I've tried to eliminate 'ufc offsets' in the code - it's a recipe for
> >>>>> disaster with renumbering.
> >>>>>
> >>>>> What we want to to know whether we have stand-alone function space or
> >>>>> a view into another space. This can't be checked at the moment. We
> >>>>> should be able to add a member function to FunctionaSpace to do
> >>>>> this, but I need to think about whether or not it will be robust. We
> >>>>> could check if a space is stand-alone or a view in the constructor
> >>>>> of a Function.
> >>>>>
> >>>>> It would be nice to mirror the design that we have for Function,
> >>>>> where we have deep and shallow copies. The idea would then be that
> >>>>> only deep copies could be used to create Functions.
> >>>>>
> >>>>> Garth
> >>>>
> >>>> We had a discussion a while back about adding a member function that
> >>>> could be used to compare the size (as in the number of dofs) of the
> >>>> subspace (sub dofmap) with the size of the bigger space. Maybe that
> >>>> could be a solution.
> >>>
> >>> That sounds usefull. But I cannot see how that can be used to check
> >>> wether a FunctionSpace is a subspace (with offset) of another one.
> >>
> >> It could if the subspace returned another value (the number of dofs)
> >> which was smaller than the whole space.
> >
> > But then the subspace need to know its parent FunctionSpace, which it afaik
> > does not know.
>
> Correct.
>
> > But if we make a SubSpace a view of another FunctionSpace it is
> > natural to add notion of which space its parent is.
> >
>
> That is tricky - the parent could go out of scope.
>
> Follow the Function design, I suggest that we require a choice at
> creation whether we want a view (SubSpace) or a fully-fledged function
> space (FunctionSpace). We can then throw an error when attempting to
> create a Function using a SubSpace.

Agree. We should avoid too fancy designs, remember the problems with
the old Function class that tried to be two things at once.

--
Anders

-- 
Problem using sub-spaces
https://bugs.launchpad.net/bugs/681415
You received this bug notification because you are a member of DOLFIN
Team, which is subscribed to DOLFIN.

Status in DOLFIN: Confirmed

Bug description:
Hi,

I have trouble using MixedFunctionSpace. Basically I want to project a function to a subspace of a MixedFunctionSpace:

   from dolfin import *
   mesh=UnitInterval(10)
   U = FunctionSpace(mesh, "CG", 1)
   V = FunctionSpace(mesh, "CG", 1)
   W=U*V
   f1=Function(W.sub(0))
   f1.vector()[:]=1.0
   #f2=project(f1, V)   # This works!
   f2=project(f1, W.sub(1))  # This doesn't!


The output of that script is attached to the end of this question. Now, the confusing thing is that the projection works fine if one projects to V instead of W.sub(1), although they are mathematically the same spaces. 

I hope I haven't missed a similar question already reported. 

Best wishes,

Simon


[0]PETSC ERROR: ------------------------------------------------------------------------
[0]PETSC ERROR: Caught signal number 11 SEGV: Segmentation Violation, probably memory access out of range
[0]PETSC ERROR: Try option -start_in_debugger or -on_error_attach_debugger
[0]PETSC ERROR: or see http://www.mcs.anl.gov/petsc/petsc-as/documentation/troubleshooting.html#Signal[0]PETSC ERROR: or try http://valgrind.org on GNU/linux and Apple Mac OS X to find memory corruption errors
[0]PETSC ERROR: configure using --with-debugging=yes, recompile, link, and run 
[0]PETSC ERROR: to get more information on the crash.
[0]PETSC ERROR: --------------------- Error Message ------------------------------------
[0]PETSC ERROR: Signal received!
[0]PETSC ERROR: ------------------------------------------------------------------------
[0]PETSC ERROR: Petsc Release Version 3.1.0, Patch 3, Fri Jun  4 15:34:52 CDT 2010
[0]PETSC ERROR: See docs/changes/index.html for recent updates.
[0]PETSC ERROR: See docs/faq.html for hints about trouble shooting.
[0]PETSC ERROR: See docs/index.html for manual pages.
[0]PETSC ERROR: ------------------------------------------------------------------------
[0]PETSC ERROR: Unknown Name on a linux-gnu named doodson by sf1409 Thu Nov 25 13:17:44 2010
[0]PETSC ERROR: Libraries linked from /build/buildd/petsc-3.1.dfsg/linux-gnu-c-opt/lib
[0]PETSC ERROR: Configure run at Fri Sep 10 04:57:14 2010
[0]PETSC ERROR: Configure options --with-shared --with-debugging=0 --useThreads 0 --with-clanguage=C++ --with-c-support --with-fortran-interfaces=1 --with-mpi-dir=/usr/lib/openmpi --with-mpi-shared=1 --with-blas-lib=-lblas --with-lapack-lib=-llapack --with-umfpack=1 --with-umfpack-include=/usr/include/suitesparse --with-umfpack-lib="[/usr/lib/libumfpack.so,/usr/lib/libamd.so]" --with-spooles=1 --with-spooles-include=/usr/include/spooles --with-spooles-lib=/usr/lib/libspooles.so --with-hypre=1 --with-hypre-dir=/usr --with-scotch=1 --with-scotch-include=/usr/include/scotch --with-scotch-lib=/usr/lib/libscotch.so --with-hdf5=1 --with-hdf5-dir=/usr
[0]PETSC ERROR: ------------------------------------------------------------------------
[0]PETSC ERROR: User provided function() line 0 in unknown directory unknown file
--------------------------------------------------------------------------
MPI_ABORT was invoked on rank 0 in communicator MPI_COMM_WORLD 
with errorcode 59.

NOTE: invoking MPI_ABORT causes Open MPI to kill all MPI processes.
You may or may not see output from other processes, depending on
exactly when Open MPI kills them.
--------------------------------------------------------------------------










References