# Result of Steady state dynamic linear elastic isnt correct

Hello,

I am new to FEnics, I want to solve the steady state dynamic linear elastic model in solid. my equation is function of frequency and the strong form is:

Divergence(

BC: Stress(vec(x),w) n(x)=T(x,w)

u(x,w)=U0

the physical problem is a plate with dimension of 1*1*0.1 with a harmonic load on the all the top and also the bottom of the plate is fixed.

I solved this as a test with the fixed frequency (w=200)

The weak form of the this equation is :

inner(sigma(u), sym(grad(v)))*dx -po*w*w *inner(

I compared the results with the model with converged mesh from Abaqus, I almost checked everything in the code the order of the displacement is mostly correct but the signs of the displacements are wrong and also the error is large between the solution of the Abaqus and this code.

I checked the variation form so I cant find any error in it. So the only thing that I suspected to be wrong is the neumann boundary condition but I cant find any bug in it.

I would be very thankful if someone can help me to figure out what is wrong in my code that I get different answers from the code.

Here is the code::

from dolfin import *

import pickle

import numpy

import csv

po=2700

w=200

Magnitude=-100

mesh = BoxMesh(0.0,0.0,0.0 , 1.0,1.0,0.1 , 10,10,10)

V = VectorFunctionS

# Defining Domain

class Bottom(SubDomain):

def inside(self, x, on_boundary):

return near(x[1], 0.0)

class Top(SubDomain):

def inside(self, x, on_boundary):

return near(x[1], 1.0)

# Initialize sub-domain instances

top = Top()

bottom = Bottom()

# Initialize mesh function for boundary domains

boundaries = FacetFunction(

boundaries.

top.mark(

bottom.

bc = DirichletBC(V, (0.0, 0.0, 0.0), boundaries,4)

# Define new measures associated with the interior domains and exterior boundaries

ds = Measure(

# Define trial and test functions

u = TrialFunction(V)

v = TestFunction(V)

f = Expression(("0.0", "scale","0.0"), w0 = w, scale = Magnitude )

# Elasticity parameters

E, nu = 69000000000, 0.3

mu = E / (2.0*(1.0 + nu))

lmbda = E*nu / ((1.0 + nu)*(1.0 - 2.0*nu))

def sigma(u):

return 2.0*mu*sym(grad(u)) + lmbda*tr(

a = inner(sigma(u), sym(grad(v)))*dx -po*w*w *inner(v,u)*dx

L = inner(f,v)*ds(2)

# Compute solution

u = Function(V)

solve(a == L, u, bc)

# Save solution in VTK format

file = File("ElasticSo

file << u

plot(u, interactive=True)

I already wrote the weak form but since I couldn't find any example or demo similar to what I want to do . does any one can help me how should i deal with this problem and solve the problem in a range of frequency or give me some tips about this problem.

Thanks

## Question information

- Language:
- English Edit question

- Status:
- Answered

- For:
- DOLFIN Edit question

- Assignee:
- No assignee Edit question

- Last query:

- Last reply:

## Can you help with this problem?

Provide an answer of your own, or ask Bahram for more information if necessary.