I have a 6x6 matrix (A) that is a function of k and H with complex symbolic entries that contain functions. I want to reduce the computational time since I will be evaluating the matrix many times for an integral that I need to compute later down the track. The code is shown below.
Clear[ω, σ, γ, L, Ms, A, α, ωM, ωH, k, \[ScriptCapitalK], d, s, dp, μ0, H, j, q]
ω = 2*Pi*19*10^9;
σ = 4.5*10^6*4;
γ = 2*Pi*2.867*10^6;
L = 30/(10^7*10^2);
Ms = 17900;
A = 1.3/10^10;
α = (A/Ms^2)*8*Pi;
ωM = γ*Ms;
ωH[H_] := γ*H
μ0 = (4*Pi)/10^7;
d = 0.0016;
\[ScriptCapitalK][k_] = Sqrt[k^2 + ω*σ*μ0*I];
q[i_, k_, H_] :=
Root[(-I)*k*μ0*σ*ω*ωM^2 -
I*k^3*α*μ0*σ*ω*ωM^2 +
k^2*ω*ωM*#1 +
I*μ0*σ*ω^2*ωM*#1 - k*ω^2*#1^2 +
k^3*α*ωM^2*#1^2 +
k^5*α^2*ωM^2*#1^2 +
I*k*α*μ0*σ*ω*ωM^2*#1^2 -
k*α*ωM^2*#1^4 - 2*k^3*α^2*ωM^2*#1^4 +
k*α^2*ωM^2*#1^6 +
k*ω^2*\[ScriptCapitalK][k]^2 -
k^3*α*ωM^2*\[ScriptCapitalK][k]^2 -
k^5*α^2*ωM^2*\[ScriptCapitalK][
k]^2 - ω*ωM*#1*\[ScriptCapitalK][k]^2 +
k*α*ωM^2*#1^2*\[ScriptCapitalK][k]^2 +
2*k^3*α^2*ωM^2*#1^2*\[ScriptCapitalK][k]^2 -
k*α^2*ωM^2*#1^4*\[ScriptCapitalK][k]^2 -
I*k*μ0*σ*ω*ωM*ωH[H] +
k*ωM*#1^2*ωH[H] +
2*k^3*α*ωM*#1^2*ωH[H] -
2*k*α*ωM*#1^4*ωH[H] -
k*ωM*\[ScriptCapitalK][k]^2*ωH[H] -
2*k^3*α*ωM*\[ScriptCapitalK][k]^2*ωH[H] +
2*k*α*ωM*#1^2*\[ScriptCapitalK][k]^2*ωH[
H] + k*#1^2*ωH[H]^2 -
k*\[ScriptCapitalK][k]^2*ωH[H]^2 & , i];
mx[i_, k_, H_, y_] := M[i]*E^(y*q[i, k, H])
my[i_, k_, H_,
y_] := -(I k mx[i, k, H,
y] (k^2 α ωM q[i, k, H]^2 - α ωM q[
i, k, H]^4 - ωM \[ScriptCapitalK][k]^2 -
k^2 α ωM \[ScriptCapitalK][
k]^2 + α ωM q[i, k, H]^2 \[ScriptCapitalK][
k]^2 + q[i, k, H]^2 ωH[H] - \[ScriptCapitalK][
k]^2 ωH[H]))/(-I μ0 σ ω ωM q[
i, k, H] + k ω q[i, k, H]^2 -
k ω \[ScriptCapitalK][k]^2 + ωM q[i, k,
H] \[ScriptCapitalK][k]^2)
hx[i_, k_, H_,
y_] := (μ0 σ ω my[i, k, H, y] q[i, k,
H] + (k mx[i, k, H, y] +
I my[i, k, H, y] q[i, k, H]) \[ScriptCapitalK][k]^2)/(
k (q[i, k, H]^2 - \[ScriptCapitalK][k]^2))
hy[i_, k_, H_, y_] := (
I (k hx[i, k, H, y] + k mx[i, k, H, y]) - my[i, k, H, y] q[i, k, H])/
q[i, k, H]
s = 0;
BoundaryConditionFree[k_, H_, y_] :=
Evaluate[{Sum[
hy[i, k, H, y] +
my[i, k, H, y] + (I*hx[i, k, H, y]*k)/Abs[k], {i, 1, 6}]}]
BoundaryConditionDielectric[k_, H_, y_] :=
Evaluate[{Sum[(hy[i, k, H, y] + my[i, k, H, y])*
Coth[Abs[k]*(d + s)] - hx[i, k, H, y]*((Abs[k]*I)/k), {i, 1,
6}]}];
j = 1;
dp = 0;
mxExchangeBoundaryCondition[k_, H_] :=
D[Sum[mx[i, k, H, y], {i, 1, 6}], y] +
dp*Sum[mx[i, k, H, y], {i, 1, 6}];
myExchangeBoundaryCondition[k_, H_] :=
D[Sum[my[i, k, H, y], {i, 1, 6}], y] +
dp*Sum[my[i, k, H, y], {i, 1, 6}];
Needs["CCompilerDriver`"]
CCompilers[]
\[ScriptCapitalA] = Compile[{{k, _Real}, {H, _Real}},
Evaluate[
Normal[CoefficientArrays[
Join[{mxExchangeBoundaryCondition[k, H] /. y -> 0},
{mxExchangeBoundaryCondition[k, H] /.
y -> L}, {myExchangeBoundaryCondition[k, H] /. y -> 0},
{myExchangeBoundaryCondition[k, H] /. y -> L},
BoundaryConditionDielectric[k, H, 0],
BoundaryConditionFree[k, H, L]],
{M[1], M[2], M[3], M[4], M[5], M[6]}][[2]]]],
Parallelization -> True, CompilationTarget -> "C"];
I've looked after the usual suspects i.e. ensuring machine precision numbers are used, functions are compiled in C and parallelization is enabled. However it still takes 0.312s to compute A.
MATLAB does this calculation almost instantly. I don't want to use it down the track though as I will need to make symbolic adjustments which will be time consuming.
M$\endgroup$