- Thank you received: 0
Time dependent source term
6 years 5 months ago #561
by bakerk
Time dependent source term was created by bakerk
Hi,
I am new to coding and am working on using FEM to solve the wave equation with a fractional damping term and am having issues with defining a time dependent source term. This is the code I have for this based on a tutorial,
time = 0
tend=5
t = Parameter(0.0)
gausspt = (1 - t + ((math.pi)**2)*(1-t+(1/2)*t**2 - (1/6)*t**3) - (t**(gamma))/(math.gamma(1-gamma)) + (t**(1-gamma))/(math.gamma(2-gamma)) -(t**(2-gamma))/(math.gamma(3-gamma)))*sin(math.pi * x) *sin(math.pi *y)
Draw(gausspt,mesh,"ft")
while time < tend-dt/2.0:
time += dt
t.Set(time)
ft = LinearForm(fes)
ft += SymbolicLFI(gausspt*v)
(gamma is a predefined constant)
When I run this it produces an error 'incompatible function arguments'. What do you suggest doing about this?
Best wishes,
Katie
I am new to coding and am working on using FEM to solve the wave equation with a fractional damping term and am having issues with defining a time dependent source term. This is the code I have for this based on a tutorial,
time = 0
tend=5
t = Parameter(0.0)
gausspt = (1 - t + ((math.pi)**2)*(1-t+(1/2)*t**2 - (1/6)*t**3) - (t**(gamma))/(math.gamma(1-gamma)) + (t**(1-gamma))/(math.gamma(2-gamma)) -(t**(2-gamma))/(math.gamma(3-gamma)))*sin(math.pi * x) *sin(math.pi *y)
Draw(gausspt,mesh,"ft")
while time < tend-dt/2.0:
time += dt
t.Set(time)
ft = LinearForm(fes)
ft += SymbolicLFI(gausspt*v)
(gamma is a predefined constant)
When I run this it produces an error 'incompatible function arguments'. What do you suggest doing about this?
Best wishes,
Katie
- christopher
- Offline
- Administrator
Less
More
- Thank you received: 101
6 years 5 months ago - 6 years 5 months ago #565
by christopher
Replied by christopher on topic Time dependent source term
There was a bug in the pow function, that affected your code, which is fixed in v6.2.1805, if you update it should work. (if you build from source, you need to check out after this fix:
github.com/NGSolve/ngsolve/commit/764a27...e58ace6086563a2a044e
)
Note as well, that the linearform is built outside the loop and only reassembled inside, this saves the time to allocate/deallocate the vectors. Your version should work as well, but is less efficient.
Best
Christopher
Note as well, that the linearform is built outside the loop and only reassembled inside, this saves the time to allocate/deallocate the vectors. Your version should work as well, but is less efficient.
Code:
from ngsolve import *
import math
from netgen.geom2d import unit_square
from time import sleep
gamma = 0.1
mesh = Mesh(unit_square.GenerateMesh(maxh=0.2))
fes = H1(mesh,order=2)
v = fes.TestFunction()
t = Parameter(0)
gausspt = (1 - t + ((math.pi)**2)*(1-t+(1/2)*t**2 - (1/6)*t**3)
-(t**(gamma))/(math.gamma(1-gamma)) + (t**(1-gamma))/(math.gamma(2-gamma))
-(t**(2-gamma))/(math.gamma(3-gamma)))*sin(math.pi * x) *sin(math.pi *y)
Draw(gausspt,mesh,"gausspt")
time = 0
tend = 5
dt = 1e-2
ft = LinearForm(fes)
ft += SymbolicLFI(gausspt*v)
while time < tend-dt/2:
time +=dt
t.Set(time)
Redraw()
sleep(0.1)
ft.Assemble()
Best
Christopher
Last edit: 6 years 5 months ago by christopher.
Time to create page: 0.127 seconds