Hi,
I got a problem concerning the implementation of complex number and how to use the workaround by building a coupled equations system. In fact, when I use basic operations to build my matrices (stiffness and mass matrix) like with dot product or grad, I can't find a method that actually works for every type of system. In the following answer, found in the Q&A, the dot product is define this way:
r_r=dot( u_r, v_r ) + dot( u_i, v_i )
r_i=dot( u_i, v_r ) - dot( u_r, v_i )
https://fenicsproject.org/qa/9136/complex-inner-product
And my dot product is slightly different but it actually gives the good result(I compared with Comsol). For the following example , the system to solve is :
$$ [K-\frac{\omega^2}{c_f^2} M +j\frac{\omega}{c_f*\rho}C] U= d$$
That I transform into a system where the j 'complex' is placed into the C matrix forming an iC matrix and I change the sign during the building of the iC matrix(ic=r_r*ds2-r_i*ds2 ). The system is now:
$$ [K-\frac{\omega^2}{c_f^2} M +\frac{\omega}{c_f*\rho}C_i] U= d $$
My problem is, for other type of linear system I don't observe correct results when I keep that philosophy of putting the 'j' in the matrices. Is there an example somewhere of the coupled system as a workaround for complex number?
sourcevalue =1.
omg=2*pi*f
cf = 343.8
u_r, u_i= TrialFunction(Vcomplex)
v_r, v_i = TestFunction(Vcomplex)
k= + (inner(grad(u_r), grad(v_r)) - inner(grad(u_i), grad(v_i))) * dx \
+ (inner(grad(u_r), grad(v_i)) + inner(grad(u_i), grad(v_r))) * dx
# for the stiffness matrix
r_r=dot( u_r, v_r ) - dot( u_i, v_i )
r_i=dot( u_i, v_r ) + dot( u_r, v_i )
m=r_r*dx +r_i*dx # for the mass matrix
c_i=r_r*ds2-r_i*ds2 #for the damping matrix
K=Matrix()
M=Matrix()
C_i=Matrix()
assemble(k, tensor=K)
assemble(m, tensor=M)
assemble(c_i, tensor=C_i)
L = ((v_r*sourcevalue)*ds1 + (v_i*sourcevalue)*ds1)
d=assemble(L)
Atot=K+(-((pow(omg,2))/(pow(cf,2)))*M)+(omg/(cf*rhof)*C_i)
u= Function(Vcomplex)
U=u.vector()
solve(Atot, U ,d)
Thanks for any type of help.
