View source

<nowiki>
<html>
##################################################################
##################################################################
#setting the domain size for the problem to be solved
domain_size = 2

##################################################################
##################################################################
## ATTENTION: here the order is important

#including kratos path
kratos_libs_path = '../../../../libs' ##kratos_root/libs
#kratos_libs_path = 'C:/kratosR1/libs' ##kratos_root/libs
kratos_applications_path = '../../../../applications/' ##kratos_root/applications
import sys
sys.path.append(kratos_libs_path)
sys.path.append(kratos_applications_path)

#importing Kratos main library
from Kratos import *
kernel = Kernel()   #defining kernel

#importing applications
import applications_interface
applications_interface.Import_IncompressibleFluidApplication = True
applications_interface.ImportApplications(kernel, kratos_applications_path)

## from now on the order is not anymore crucial
##################################################################
##################################################################

from KratosR1IncompressibleFluidApplication import *


#defining a model part
model_part = ModelPart("FluidPart");

##importing the solver files and adding the variables
import incompressible_fluid_solver
incompressible_fluid_solver.AddVariables(model_part)

#adding of Variables to Model Part should be here when the "very fix container will be ready"

#reading a model
gid_io = GidIO("cavity2d",GiDPostMode.GiD_PostBinary)
##gid_io.ReadMesh(model_part.GetMesh())
gid_io.ReadModelPart(model_part)
gid_io.WriteMesh((model_part).GetMesh(),domain_size,GiDPostMode.GiD_PostBinary);
print model_part

#the buffer size should be set up here after the mesh is read for the first time
model_part.SetBufferSize(3)

##add Degrees of Freedom to all of the nodes
incompressible_fluid_solver.AddDofs(model_part)



#creating a fluid solver object
fluid_solver = incompressible_fluid_solver.IncompressibleFluidSolver(model_part,domain_size)
fluid_solver.laplacian_form = 2;
fluid_solver.predictor_corrector = True
fluid_solver.vel_toll = 1e-3
fluid_solver.time_order = 2
fluid_solver.echo_level = 0



#pILUPrecond = ILU0Preconditioner() 
#fluid_solver.pressure_linear_solver =  BICGSTABSolver(1e-9, 5000,pILUPrecond)
pDiagPrecond = DiagonalPreconditioner()
fluid_solver.velocity_linear_solver =  BICGSTABSolver(1e-9, 5000,pDiagPrecond)
fluid_solver.pressure_linear_solver =  BICGSTABSolver(1e-9, 5000,pDiagPrecond)
##fluid_solver.pressure_linear_solver = SkylineLUFactorizationSolver();
fluid_solver.Initialize()

#settings to be changed
Re = 100.0
nsteps = 200
output_step = 1

Dt = 1000.0/Re
if(Dt > 0.1):
    Dt = 0.1
out = 0
for node in model_part.Nodes:
    node.SetSolutionStepValue(VISCOSITY,0,1.0/Re)

for step in range(1,nsteps):
    print "line49"

    time = Dt*step
    print time
    model_part.CloneTimeStep(time)

    print "qui"

    print time
    #print model_part.ProcessInfo()[TIME]

    #solving the fluid problem
    if(step > 3):
        fluid_solver.Solve()

    print "li"


    #print the results
    if(out == output_step):
        gid_io.WriteNodalResults(PRESSURE,model_part.Nodes,time,0)
        gid_io.WriteNodalResults(VELOCITY,model_part.Nodes,time,0)
        gid_io.WriteNodalResults(PRESS_PROJ,model_part.Nodes,time,0)
        gid_io.WriteNodalResults(CONV_PROJ,model_part.Nodes,time,0)
        gid_io.WriteNodalResults(NODAL_AREA,model_part.Nodes,time,0)
        gid_io.WriteNodalResults(VISCOSITY,model_part.Nodes,time,0)
        out = 0
    out = out + 1

node = model_part.Nodes[1]

print node
</html>
</nowiki>