Problem with .Assemble()

Main Forum Old-Forum - do not post here
1

Hello everyone,

I am trying to implement a k-epsilon model for a turbulent flow in a three dimensional geometry.
When I try to run code, the program runs up to an assemble statement and then just exits without giving me an error message. Does someone have any insights on why the program just breaks doing an assemble statement.
Here is my code.

from ngsolve import *
from netgen import NgOCC
from math import *
import sys


geo = NgOCC.LoadOCCGeometry("Normale_Geometrie_3.step")
mesh= Mesh(geo.GenerateMesh())

dirichletstring = "bc_0|bc_1|bc_2"
for i in range(4,103):
    dirichletstring = dirichletstring + "|" + "bc_" + "%d" % (i)
    
V = VectorH1(mesh,order=3, dirichlet=dirichletstring)
Q = H1(mesh,order=2)
Q_alt = H1(mesh, order=2, dirichlet = dirichletstring)
N = NumberSpace(mesh)

X = FESpace([V,Q])

#constants

density = 0.006984526006623
sigma_k = 1
sigma_e = 1.3
C_e1 = 1.44
C_e2 = 1.92
C_mu = 0.09

#constant viscosity

nu = 35.1524056250947

#relaxation parameters

omega_k = 0.4
omega_T = 0.4
omega_e = 0.4

#trial and test functions

u,p = X.TrialFunction()
v,q = X.TestFunction()
u_star = V.TrialFunction()
phi = Q.TrialFunction()
l= V.TrialFunction()
k_star = Q_alt.TrialFunction()
e_star = Q_alt.TrialFunction()
e = Q_alt.TrialFunction() 
v_new = Q_alt.TestFunction()

#gridfunctions

u_star_gfu = GridFunction(V)
phi_gfu = GridFunction(Q)
l_gfu = GridFunction(V)
k_star_gfu = GridFunction(Q_alt)
k_gfu_last = GridFunction(Q_alt)
k_gfu_next = GridFunction(Q_alt)
k_iter_gfu_last = GridFunction(Q_alt)
k_iter_gfu_next = GridFunction(Q_alt)
e_star_gfu = GridFunction(Q_alt)
e_gfu_last = GridFunction(Q_alt)
e_gfu_next = GridFunction(Q_alt)
e_iter_gfu_last = GridFunction(Q_alt)
e_iter_gfu_next = GridFunction(Q_alt)

p_k_gfu = GridFunction(Q_alt)

f_2_gfu = GridFunction(Q_alt)
f_mu_gfu = GridFunction(Q_alt)
Re_T_gfu = GridFunction(Q_alt)

gfu_last = GridFunction(X)
gfu_next = GridFunction(X)

nu_T_star_gfu = GridFunction(N)
nu_T_gfu_last = GridFunction(N)
nu_T_gfu_next = GridFunction(N)
nu_T_iter_gfu_last = GridFunction(N)
nu_T_iter_gfu_next = GridFunction(N)

#boundary conditions

u_inlet = CoefficientFunction( (0,0,5))
gfu_last.components[0].Set(u_inlet, definedon=mesh.Boundaries("bc_103"))

k_inlet = 0.01 * 25   
e_inlet = C_mu * sqrt(k_inlet) * k_inlet / 1100 


k_star_gfu.Set(CoefficientFunction(k_inlet), definedon = mesh.Boundaries("bc_103"))
k_gfu_last.Set(CoefficientFunction(k_inlet), definedon = mesh.Boundaries("bc_103"))
k_gfu_next.Set(CoefficientFunction(k_inlet), definedon = mesh.Boundaries("bc_103")) 
k_iter_gfu_last.Set(CoefficientFunction(k_inlet), definedon = mesh.Boundaries("bc_103"))
k_iter_gfu_next.Set(CoefficientFunction(k_inlet), definedon = mesh.Boundaries("bc_103"))
e_star_gfu.Set(CoefficientFunction(e_inlet), definedon = mesh.Boundaries("bc_103"))
e_gfu_last.Set(CoefficientFunction(e_inlet), definedon = mesh.Boundaries("bc_103"))
e_gfu_next.Set(CoefficientFunction(e_inlet), definedon = mesh.Boundaries("bc_103"))
e_iter_gfu_last.Set(CoefficientFunction(e_inlet), definedon = mesh.Boundaries("bc_103"))
e_iter_gfu_next.Set(CoefficientFunction(e_inlet), definedon = mesh.Boundaries("bc_103"))


#initial conditions 

mixing_length = 1100

k_0 = (nu/mixing_length)* (nu/mixing_length)
e_0 = C_mu * sqrt(k_0) * k_0 / mixing_length
Re_T_0 = k_0 * k_0 / (nu * e_0)
f_mu_0 = exp(-3.4/((1+Re_T_0/50)*(1+Re_T_0/50)))


k_gfu_last.Set(CoefficientFunction(k_0))
k_iter_gfu_last.Set(CoefficientFunction(k_0))
k_iter_gfu_next.Set(CoefficientFunction(k_0))

e_gfu_last.Set(CoefficientFunction(k_0))
e_iter_gfu_last.Set(CoefficientFunction(k_0))
e_iter_gfu_next.Set(CoefficientFunction(k_0))

nu_T_gfu_last.Set(CoefficientFunction(C_mu * f_mu_0 * (k_0 * k_0)/e_gfu_last))



 


#time specification

time= 0.0
dt= 0.1
timeend = 0.5

Draw(gfu_last.components[0],mesh,"u")


while time < timeend:
    print("time = ", time)
    #equation for u_star with beta=0

    u_star_transient_i = InnerProduct(u_star,v)
    
    a_u_star = BilinearForm(V)
    a_u_star += (u_star_transient_i +   InnerProduct(grad(u_star),grad(v)))*dx
    a_u_star.Assemble()

At this assemble statement the program just exits and I don’t know why.

Thanks for your help.

Leon

2

Issue resolved.