project conditional statement to print values makes code slow

11 months ago by
Hi, I am using conditional statement and to print the values I am projecting. It makes my code very slow. Projection is just to veryfyanswer.
from dolfin import *
import matplotlib.pyplot as plt
E = 29e3
K = 500.0
sy = 36.0
Area = 0.1
print 'properties:', E,K,sy, Area

mesh = IntervalMesh(60,0,60)
#plot(mesh, interactive = True)

V  = FunctionSpace(mesh,'CG',1)
P1 = FunctionSpace(mesh,'DG',0)
u = TrialFunction(V)
v = TestFunction(V)

F = E*inner(grad(u),grad(v))*dx

a = lhs(F)
L = rhs(F)

u = Function(V)

class right(SubDomain):
    def inside(self,x,on_boundary):
        tol = 1e-10
        return abs(x[0]-60.0) < tol and on_boundary
class left(SubDomain):
    def inside(self,x,on_boundary):
        tol = 1e-10
        return abs(x[0]-0.0) < tol and on_boundary
Right = right()
Left = left()
uml = Expression('t', t = 0.0)
bcl= DirichletBC(V, Constant(0.0), Left)
bcr = DirichletBC(V, uml, Right)

bc = [bcl, bcr]
t = 0.0
dt = 0.005
num_steps = (0.5-0.0)/dt
num_steps = 100
delta = 0.00005
pStrain = 0.0
alphan = 0.0
pStrainN = 0.0

def stre(Gamma,E,Tstress):
    return (1.0-(Gamma*E)/abs(Tstress))*Tstress
def E_val(K,E):
    return E*K/(E+K)
def alpha_val(alphan,Gamma):
    return alphan+Gamma
def pStrain_val(pStrainN,Gamma,Tstress ):
    return pStrainN+ Gamma*sign(Tstress)

for n in range(num_steps):
    uml.t = t
    solve(a==L, u, bc)
    eps = u.dx(0)
    #print u.vector().norm('l2')
    ep = project(eps,P1)
    Tstrain = eps + delta
    Tstress = E*(Tstrain - pStrain)
    pStrain = pStrainN
    fun = abs(Tstress) -(sy + K*alphan)
    #fu = project(fun,P1)    
    stress = Tstress    
    st = project(E*(Tstrain),P1)
    right = Point(60.0)
    print 'str:',ep(right), 'Tstress', st(right)
    Gamma = fun/(E+K)
    stress = conditional(le(fun,0.0),Tstress,stre(Gamma,E,Tstress))
    C_ep = conditional(le(fun,0.0),E, E_val(K,E) )
    pStrain = conditional(le(fun,0.0),pStrain,pStrain_val(pStrainN,Gamma,Tstress ) )
    alphaT = conditional(le(fun,0.0),alphan,alpha_val(alphan,Gamma))
    alphan = alphaT
    pStraiN = pStrain
    if n>=20:
        stt = project(stress,P1)
        print 'stress:', stt(right)​

Any suggestions.

Thanks for the help
Community: FEniCS Project
You could start by reducing your code to a minimal test case. For example, I see a time loop, and many print statements. How can you ask to us to spend our time on _your_ problem if _you_ are not willing to spend time on it? Since projection is to verify the result, you could likely drop all the  code needed to compute Tstress and still show your problem.
After this, you could also repost your question with some timings, with and without the conditional statement, so that who reads you can understand what you mean by "very slow", and compare the timing with and without projection.
written 11 months ago by Marco Morandini  

1 Answer

11 months ago by
Did you try to use an iterative solver the projection? I think FEniCS recently changed the default behaviour of projection from CG+ILU to a direct solver.
Please login to add an answer/comment or follow this question.

Similar posts:
Search »