No commit message

No commit message

Solver/TESTCASES/Advection1D.lua

+
+--[[ 
+     Function for initial conditions
+--]]
+
+function initial_condition( xyz , f )
+  for i=0,xyz:size1()-1 do
+    x = xyz:get(i,0)
+    y = xyz:get(i,1)
+    z = xyz:get(i,2)
+    if (x>-0.3 and x<0.3) then
+     f:set (i, 0, 1)
+   else
+     f:set (i, 0, 0)
+    end	
+  end
+end
+
+--[[ 
+     Example of a lua program driving the DG code
+--]]
+
 model = GModel  ()
-model:load ('edge.geo')
--- model:mesh (1)
--- model:save ('edge.msh')
 model:load ('edge.msh')
-dg = dgSystemOfEquations (model)
-dg:setOrder(1)
+order=1
+dimension=1
 
 -- conservation law
 -- advection speed
@@ -17,9 +36,6 @@ nu=fullMatrix(1,1);
 nu:set(0,0,0)
 
 law = dgConservationLawAdvectionDiffusion(functionConstant(v):getName(), '') 
---FunctionConstant(nu):getName())
-
-dg:setConservationLaw(law)
 
 -- boundary condition
 outside=fullMatrix(1,1)
@@ -28,35 +44,31 @@ bndcondition=law:newOutsideValueBoundary(functionConstant(outside):getName())
 law:addBoundaryCondition('Left',bndcondition)
 law:addBoundaryCondition('Right',bndcondition)
 
-dg:setup()
+groups = dgGroupCollection(model, dimension, order)
+groups:buildGroupsOfInterfaces()
 
--- initial condition
-function initial_condition( xyz , f )
-  for i=0,xyz:size1()-1 do
-    x = xyz:get(i,0)
-    y = xyz:get(i,1)
-    z = xyz:get(i,2)
-    if (x>-0.3 and x<0.3) then
-     f:set (i, 0, 1)
-   else
-     f:set (i, 0, 0)
-    end	
-  end
-end
-dg:L2Projection(functionLua(1,'initial_condition',{'XYZ'}):getName())
-print'***exporting init solution ***'
+-- build Runge Kutta and limiter
+rk=dgRungeKutta()
+limiter = dgSlopeLimiter(law)
+rk:setLimiter(limiter) 
 
-dg:exportSolution('output/Adv1D_unlimited')
-dg:limitSolution()
-dg:exportSolution('output/Adv1D-00000')
+-- build solution vector
+FS = functionLua(1, 'initial_condition', {'XYZ'}):getName()
+solution = dgDofContainer(groups, law:getNbFields())
+solution:L2Projection(FS)
 
+solution:exportMsh('output/init')
+limiter:apply(solution)
+solution:exportMsh('output/init_limit')
 
--- main loop
+print'*** solve ***'
+local x = os.clock()
 n = 5
+dt = 0.03
 for i=1,100*n do
-  norm = dg:RK44_limiter(0.03)
+  norm = rk:iterate44(law,dt,solution)
   if (i % n == 0) then 
-    print('iter',i,norm)
-    dg:exportSolution(string.format("output/Adv1D-%05d", i)) 
+    print('|ITER|',i,'|NORM|',norm,'|DT|',dt,'|CPU|',os.clock() - x)
+    solution:exportMsh(string.format('output/solution-%06d',i))
   end
 end