Merged develop, solved conflicts.

src/bvpy/domains/abstract.py

@@ -104,20 +104,20 @@ class AbstractDomain(ABC):
     bdata : :class:`MeshFunctionSizet<dolfin.cpp.mesh.MeshFunctionSizet>`
         Description of attribute `bdata`.
     dx : :class:`Measure<ufl.measure.Measure>`
-        Description of attribute `dx`.
-    ds : type
-        Description of attribute `ds`.
-    mesh : type
-        Description of attribute `mesh`.
+        An integral measure within the domain of interest.
+    ds : :class:`Measure<ufl.measure.Measure>`
+        An integral measure on the border of the domain of interest.
+    mesh : :class:`Mesh<dolfin.cpp.mesh.Mesh>`
+        A mesh representation (topology and geometry) of the considered domain.
     model : :class:`Model<gmsh.model>`
-        The type of mesh model used ????
+        A purely geometrical description of the domain.
     model_name : str
-        Name of the model.
+        Name of the geometrical model of the domain.
 
     Methods
     -------
     _mesh_conversion
-        Converts a mesh into ????
+        Converts a temporary gmsh mesh into a fenics (dolphin) one.
     _set_factory
         Sets the factory used to generate the domain elements?
     geometry
@@ -220,8 +220,10 @@ class AbstractDomain(ABC):
 
         Yields
         -------
-        ???
-            Description of returned object.
+        :class:`occ<gmsh.model.occ>`
+        or
+        :class:`occ<gmsh.model.geo>`
+            A gmsh model factory to generate mesh from geometry.
 
         """
         if factory == 'occ':
@@ -247,7 +249,7 @@ class AbstractDomain(ABC):
 
     @property
     def model_name(self):
-        """Names the model.
+        """Name of the geometrical model of the domain.
 
         Returns
         -------
@@ -259,18 +261,18 @@ class AbstractDomain(ABC):
 
     @property
     def model(self):
-        """The gmsh model used within the domain.
+        """The (gmsh) geometrical model used to define the domain.
 
         Returns
         -------
         :class:`Model<gmsh.model>`
-            The gmsh model used to tile the domain.
+            A purely geometrical description of the domain.
 
         """
         return AbstractDomain._model
 
     def synchronize(self):
-        """synchronizes the element factory ????
+        """synchronizes the CAD model with the gmsh model.
 
         Returns
         -------
@@ -286,10 +288,7 @@ class AbstractDomain(ABC):
 
         Parameters
         ----------
-        *args : type
-            Description of parameter `*args`.
-        **kwargs : type
-            Description of parameter `**kwargs`.
+        None
 
         Returns
         -------
@@ -303,7 +302,7 @@ class AbstractDomain(ABC):
         """
 
     def run_gui(self):
-        """???
+        """Generates a quick gmsh gui with the current model in it.
 
         Returns
         -------
@@ -351,16 +350,29 @@ class AbstractDomain(ABC):
             gmsh.option.setNumber("General.Verbosity", 0)
 
     def generate(self):
-        """Generates the data structure and writes it on disc.
+        """Generates a fenics mesh representing the domain.
 
         Yields
         ------
-        :class:`Mesh<gmsh.mesh>` ???
-            Tile mesh corresponding to the tiled domain.
-        :class:'ds<fenics.mesh.MeshFunction>' ??
-            The boundary surface incremental element.
-        :class:'dx<fenics.mesh.MeshFunction>' ??
-            The volume incremental element. ???
+        mesh : :class:`Mesh<dolfin.cpp.mesh.Mesh>`
+            A fenics mesh corresponding to the domain.
+        bdata : :class:`MeshFunctionSizet<dolfin.cpp.mesh.MeshFunction>`
+            A meshFunction that labels the borders (i.e. outermost facets)
+            of the mesh of the domain.
+        cdata : :class:`MeshFunctionSizet<dolfin.cpp.mesh.MeshFunction>`
+            A meshFunction that labels the inner elements (i.e. cells)
+            of the mesh of the domain.
+        ds : :class:`Measure<ufl.measure.Measure>`
+            An integral measure on the border of the domain of interest.
+        dx : :class:`Measure<ufl.measure.Measure>`
+            An integral measure within the domain of interest.
+
+        Notes
+        -----
+        * The method also generates a temporary gmsh mesh
+        and stores it temporary on disc.
+        * MeshFunctionSizet : `Sizet` means that the output of
+        the function is a positive integer.
 
         """
         if self._cell_type == 'line':
@@ -380,12 +392,18 @@ class AbstractDomain(ABC):
         self.dx = fe.Measure('dx', domain=self.mesh, subdomain_data=self.cdata)
 
     def _mesh_conversion(self):
-        """Short summary.
+        """Converts a temporary gmsh mesh into a fenics (dolphin) one.
 
         Returns
         -------
-        type
-            Description of returned object.
+        mesh : :class:`Mesh<dolfin.cpp.mesh.Mesh>`
+            A fenics mesh corresponding to the domain.
+        bdata : :class:`MeshFunctionSizet<dolfin.cpp.mesh.MeshFunction>`
+            A meshFunction that labels the borders (i.e. outermost facets)
+            of the mesh of the domain.
+        cdata : :class:`MeshFunctionSizet<dolfin.cpp.mesh.MeshFunction>`
+            A meshFunction that labels the inner elements (i.e. cells)
+            of the mesh of the domain.
 
         """
         path = self._dir_name + '/mesh.msh'
@@ -443,7 +461,8 @@ class AbstractDomain(ABC):
                         line_data =\
                             mesh_io.cell_data_dict["gmsh:physical"][key]
                     else:
-                        line_data = np.vstack([line_data, mesh_io.cell_data_dict["gmsh:physical"][key]])
+                        line_data = np.vstack([line_data,
+                                mesh_io.cell_data_dict["gmsh:physical"][key]])
 
             line_mesh = meshio.Mesh(points=mesh_io.points,
                                     cells=[("line", line_cells)],
@@ -460,12 +479,62 @@ class AbstractDomain(ABC):
 
 
 class BuiltInModel(object):
+    """A model based on the Gmsh `geo` factory.
+
+    Parameters
+    ----------
+    None
+
+    Attributes
+    ----------
+    factory : :class:`occ<gmsh.model.geo>`
+        A gmsh model `geo` factory to generate mesh from geometry.
+
+    """
+
     def __init__(self, *args, **kwargs):
+        """Sets the attribute factory to a specific type of gmsh factory.
+
+        Parameters
+        ----------
+        None
+
+        Returns
+        -------
+        factory : :class:`occ<gmsh.model.geo>`
+            A gmsh model `geo` factory to generate mesh from geometry.
+
+        """
         super().__init__()
         self.factory = gmsh.model.geo
 
 
 class OccModel(object):
+    """A model based on the Gmsh `Occ` factory.
+
+    Parameters
+    ----------
+    None
+
+    Attributes
+    ----------
+    factory : :class:`occ<gmsh.model.occ>`
+        A gmsh model `occ` factory to generate mesh from geometry.
+
+    """
+
     def __init__(self, *args, **kwargs):
+        """Sets the attribute factory to a specific type of gmsh factory.
+
+        Parameters
+        ----------
+        None
+
+        Returns
+        -------
+        factory : :class:`occ<gmsh.model.occ>`
+            A gmsh model `occ` factory to generate mesh from geometry.
+
+        """
         super().__init__()
         self.factory = gmsh.model.occ

src/bvpy/domains/geometry.py

@@ -22,7 +22,6 @@
 
 import fenics as fe
 import numpy as np
-from math import pi
 
 
 class Init_orientation_3D(fe.UserExpression):
@@ -30,20 +29,37 @@ class Init_orientation_3D(fe.UserExpression):
 
     Attributes
     ----------
-    mesh : :class:`dolfin.cpp.mesh.Mesh<dolfin.cpp.mesh.Mesh>`
+    mesh : :class:`Mesh<dolfin.cpp.mesh.Mesh>`
         Mesh to work with.
-    interior : list(float)
+    barycenter : numpy.array of floats
+        The (3D) position vector of the barycenter of all the mesh vertices.
+    interior : list of floats
         A 3D vector (x, y, z) to orient properly the mesh.
 
     Methods
     -------
-    eval(value, x)
-        Compute the relative position with respect to the barycenter.
-    value_shape()
-        Return the geometrical dimesion of the structure (always 3).
+    eval
+        Computes the relative position with respect to the barycenter.
+    value_shape
+        Returns the geometrical dimesion of the structure (always 3).
     """
 
     def __init__(self, mesh, interior, **kwargs):
+        """Short summary.
+
+        Parameters
+        ----------
+        mesh : :class:`Mesh<dolfin.cpp.mesh.Mesh>`
+            Mesh to work with.
+        interior : list of floats
+            A 3D vector (x, y, z) to orient properly the mesh.
+
+        Returns
+        -------
+        barycenter : numpy.array of floats
+            The (3D) position vector of the barycenter of the mesh vertices.
+
+        """
         if interior is None:
             self.barycenter = mesh.coordinates().mean(axis=0)
         else:
@@ -51,11 +67,38 @@ class Init_orientation_3D(fe.UserExpression):
         super().__init__(**kwargs)
 
     def eval(self, value, x):
+        """Computes the relative position with respect to the barycenter.
+
+        Parameters
+        ----------
+        value : numpy.array of float
+            the position we want to update (???).
+        x : list of float
+            The position vector to re-center.
+
+        Returns
+        -------
+        None
+
+        Notes
+        -----
+        This method is used by fenics and is not explicitly used within
+        the bvpy code but is nonetheless needed.
+
+        """
         value[0] = x[0] - self.barycenter[0]
         value[1] = x[1] - self.barycenter[1]
         value[2] = x[2] - self.barycenter[2]
 
     def value_shape(self):
+        """Returns the geometrical dimesion of the structure (always 3).
+
+        Returns
+        -------
+        tuple of int
+            This always equals 3.
+
+        """
         return (3,)
 
 
@@ -64,13 +107,17 @@ def normal(domain, scale=1, interior=None, degree=1, name='normal'):
 
     Parameters
     ----------
-    domain : Domain
+    domain : subclass of :class:`Domain<bvpy.domains.AbstracDomain>`
         The meshed surface on which we want to compute the normal field.
-    scale : type
-        Optional (Default : 1). The displacement intensity to consider.
-    interior : list(float)
-        Optional (Default : None). A 3D vector (x, y, z)
-        to orient properly the mesh.
+    scale : float, optional
+        The displacement intensity to consider (the default is 1).
+    interior : list(float), optional
+        A 3D vector to set the mesh orientation (the default is None).
+    degree : int, optional
+        The interpolation degree to use to compute the vector space
+        on the border (the default is 1).
+    name : str, optional
+        the label we want to attach to the normal field.
 
     Returns
     -------
@@ -92,7 +139,7 @@ def normal(domain, scale=1, interior=None, degree=1, name='normal'):
     V_scal = fe.FunctionSpace(mesh, 'CG', 1)
     vf = np.zeros([mesh.num_entities(0), 3])
 
-    mesh.init(0,2)
+    mesh.init(0, 2)
 
     for v in fe.vertices(mesh):
         f_vertex = np.zeros(3)
@@ -121,29 +168,75 @@ def normal(domain, scale=1, interior=None, degree=1, name='normal'):
     return u
 
 def get_boundary_facets(mesh):
+    """Selects the boundary elements of a given mesh.
+
+    Parameters
+    ----------
+    mesh : :class:`Mesh<dolfin.cpp.mesh.Mesh>`
+        The mesh from which we want to get the boundary elements.
+
+    Returns
+    -------
+    facets : list of :class:`facets<dolfin.cpp.mesh.facets>`
+        the list of the outermost elements of the considered mesh.
+
+    Notes
+    -----
+    The returned elements corresponds to mesh elements of topological
+    dimension (n-1) for a mesh of topological dimension n.
+
+    """
     mesh.init()
     facets = list()
 
     for fct in fe.facets(mesh):
-        if len(list(fe.cells(fct)))==1:
+        if len(list(fe.cells(fct))) == 1:
             facets.append((fct))
+
     return facets
 
 
 def get_boundary_vertex(mesh):
+    """Short summary.
+
+    Parameters
+    ----------
+    mesh : :class:`Mesh<dolfin.cpp.mesh.Mesh>`
+        The mesh from which we want to get the boundary vertices.
+
+
+    Returns
+    -------
+    list of :class:`Vertex<dofin.cpp.mesh.Vertex>`
+        The list of the outermost vertices of the considered mesh.
+
+    """
     bnd_facets = get_boundary_facets(mesh)
     bnd_vertex = list()
 
     for fct in bnd_facets:
         [bnd_vertex.append(v.index()) for v in fe.vertices(fct)]
 
-    bnd_vertex = list(set(bnd_vertex)) #Unique
+    bnd_vertex = list(set(bnd_vertex))  # Unique
     bnd_vertex = [fe.Vertex(mesh, i) for i in bnd_vertex]
 
     return bnd_vertex
 
 
 def boundary_normal(mesh):
+    """Computes the normal vector field on the boundary of a mesh.
+
+    Parameters
+    ----------
+    mesh : :class:`Mesh<dolfin.cpp.mesh.Mesh>`
+        The mesh on which the boundary normals should be computed.
+
+    Returns
+    -------
+    :class:`Function<dolfin.functions.function.Function>`
+        a list of 3D vector defined on the outermost vertices of a mesh.
+
+    """
     mesh.init()
     bnd_facets = get_boundary_facets(mesh)
     facets_index = [f.index() for f in bnd_facets]
@@ -151,9 +244,11 @@ def boundary_normal(mesh):
 
     for fct in bnd_facets:
         for vtx in fe.vertices(fct):
-            norm = np.array([f.normal().array() for f in fe.facets(vtx) if f.index() in facets_index] )
+            norm = np.array([f.normal().array()
+                             for f in fe.facets(vtx)
+                             if f.index() in facets_index])
             norm = np.sum(norm, axis=0)
-            norm = norm/np.linalg.norm(norm)
+            norm = norm / np.linalg.norm(norm)
             normals[vtx.index()] = norm
 
     V = fe.VectorFunctionSpace(mesh, 'P', 1, dim=3)

src/bvpy/ivp.py

@@ -33,22 +33,44 @@ class Time():
 
 
 class Scheme():
-    def __init__(self, vform):
+    def __init__(self, vform=None, dt=None):
         self.vform = vform
+        self.previous_solution = None
+        self.dt = dt
 
-    def implicit_euler(lhs, rhs):
-        return lhs, Zero()
+    def set_previous_solution(self, previous):
+        self.previous_solution = previous
+
+    def implicit_euler(self, u, v, sol):
+        lhs, rhs = self.vform.get_form(u, v, sol)
+
+        lhs = fe.Constant(self.dt)*lhs
+        rhs = fe.Constant(self.dt)*rhs
+
+        lhs += u*v*self.dx
+        rhs += self.previous_solution*v*self.dx
+
+        return lhs, rhs
 
     def explicit_euler(lhs, rhs):
         return lhs, rhs
 
+    @property
+    def dx(self):
+        return self.vform.dx
+
+    @property
+    def ds(self):
+        return self.vform.ds
+
 
 
 class IVP(BVP):
-    def __init__(self, domain=None, vform=None, bc=None):
+    def __init__(self, domain=None, vform=None, bc=None,):
         super().__init__(domain=domain, vform=vform, bc=bc)
 
         self.time = Time()
+        self.scheme = Scheme()
 
     def set_time(self, time):
         self.time.current = time

tutorials/bvpy_tutorial_3_poissonEquation_FEniCSlike.ipynb

@@ -223,4 +223,4 @@
  },
  "nbformat": 4,
  "nbformat_minor": 4
-}
\ No newline at end of file
+}

tutorials/bvpy_tutorial_5_linearElasticity_2DFlatCellularizedTissue.ipynb

@@ -62,6 +62,18 @@
   {
    "cell_type": "code",
    "execution_count": null,
+<<<<<<< HEAD
+=======
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "type(drchlt)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+>>>>>>> develop
    "metadata": {},
    "outputs": [],
    "source": [
@@ -75,6 +87,18 @@
   {
    "cell_type": "code",
    "execution_count": null,
+<<<<<<< HEAD
+=======
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "p_force.vector().get_local().reshape(-1,3)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+>>>>>>> develop
    "metadata": {},
    "outputs": [],
    "source": [