fixed sign error in rho form

no_c_fixed_bdry_fenics.py

@@ -19,6 +19,7 @@ class Growth(object):
         # define function space with this mixed element
         # define function space with this mixed element
         self.function_space = df.FunctionSpace(self.mesh, mixed_element)
         self.function_space = df.FunctionSpace(self.mesh, mixed_element)
 
 
+
         # dirichlet boundaries for u, v
         # dirichlet boundaries for u, v
         bc1 = df.DirichletBC(self.function_space.sub(0), df.Constant(0.0), 'on_boundary')   # u,v at boundary = 0
         bc1 = df.DirichletBC(self.function_space.sub(0), df.Constant(0.0), 'on_boundary')   # u,v at boundary = 0
         bc2 = df.DirichletBC(self.function_space.sub(1), df.Constant(0.0), 'on_boundary')
         bc2 = df.DirichletBC(self.function_space.sub(1), df.Constant(0.0), 'on_boundary')
@@ -33,7 +34,7 @@ class Growth(object):
         return (df.inner((vel * conc).dx(0), tconc))
         return (df.inner((vel * conc).dx(0), tconc))
 
 
     def reaction_rho(self, rho, trho):
     def reaction_rho(self, rho, trho):
-        return df.inner(self.turnover_rho * (rho - self.average_rho), trho)
+        return self.turnover_rho * df.inner(rho - self.average_rho, trho)
 
 
     def  diffusion_reaction_rho(self, rho, trho):
     def  diffusion_reaction_rho(self, rho, trho):
         return (self.diffusion_rho * df.inner(rho.dx(0), trho.dx(0)) + self.reaction_rho(rho, trho))
         return (self.diffusion_rho * df.inner(rho.dx(0), trho.dx(0)) + self.reaction_rho(rho, trho))
@@ -43,7 +44,9 @@ class Growth(object):
         u0 = df.interpolate(u0, self.function_space.sub(0).collapse())
         u0 = df.interpolate(u0, self.function_space.sub(0).collapse())
         rho0 = df.interpolate(rho0, self.function_space.sub(2).collapse())
         rho0 = df.interpolate(rho0, self.function_space.sub(2).collapse())
 
 
+        # v0_function_space = df.FunctionSpace(self.mesh, 'P', 1, constrained_domain = self.pbc)
         v0_function_space = df.FunctionSpace(self.mesh, 'P', 1)
         v0_function_space = df.FunctionSpace(self.mesh, 'P', 1)
+
         v0  = df.Function(v0_function_space)
         v0  = df.Function(v0_function_space)
         tv0 = df.TestFunction(v0_function_space)
         tv0 = df.TestFunction(v0_function_space)
         v0form = (df.inner(v0,tv0)
         v0form = (df.inner(v0,tv0)
@@ -82,8 +85,8 @@ class Growth(object):
         )
         )
 
 
         rhoform = (df.inner((rho - rho0)/self.timestep, trho)
         rhoform = (df.inner((rho - rho0)/self.timestep, trho)
-                - 3/2 * self.advection(rho0, v, trho)
-                + 1/2 * self.advection(rhominus1, v, trho)
+                + 3/2 * self.advection(rho0, v, trho)
+                - 1/2 * self.advection(rhominus1, v, trho)
                 + 9/16 * self.diffusion_reaction_rho(rho, trho)
                 + 9/16 * self.diffusion_reaction_rho(rho, trho)
                 + 3/8 * self.diffusion_reaction_rho(rho0, trho)
                 + 3/8 * self.diffusion_reaction_rho(rho0, trho)
                 + 1/16 * self.diffusion_reaction_rho(rhominus1, trho)
                 + 1/16 * self.diffusion_reaction_rho(rhominus1, trho)
@@ -128,8 +131,12 @@ if __name__ == '__main__':
 
 
     g = Growth(params)
     g = Growth(params)
 
 
-    u0 = df.Expression('sin(x[0])', degree=1)  
-    rho0 = df.Expression('cos(x[0])+cos(x[0]/2)', degree=1)
+    u0 = df.Expression('sin(x[0])', degree=1)
+    u0 = df.interpolate(u0, g.function_space.sub(0).collapse())
+
+    rho0 = df.Expression('cos(x[0])+cos(x[0]/2)', degree=1) 
+    rho0 = df.interpolate(rho0, g.function_space.sub(2).collapse())
+
     (u, v, rho) = g.solve(u0, rho0)
     (u, v, rho) = g.solve(u0, rho0)
 
 
     x = g.mesh.coordinates()
     x = g.mesh.coordinates()