1
$\begingroup$

In my programme I have defined all these functions:

getA[kappa_] := 
  Table[2*Cos[(2*Pi/n)*(Abs[j - i])*kappa], {i, n}, {j, n}];
getF[csi_, a_, b_] := 
  Module[{csiInv = Inverse[csi]}, .5 Tr[csiInv.a.csiInv.b]];
getG[csi_, f_, a_] := Module[{csiInv = Inverse[csi]}, 
   csiInv.a.csiInv/2];
getE[g_, k_] := 
  Module[{kinv = Inverse[k], ktransinv2 = Transpose[kinv]},
   ktransinv2.g.kinv];
getW[k_, a_, e_] := Module[{ktrans = Transpose[k]},
   Tr[k.a.ktrans.e]];
getV[csi_, e_, e2_, k_] := 
  Module[{ktrans = Transpose[k], e2trans = Transpose[e2]},
   2*Tr[csi.ktrans.e.k.csi.ktrans.e2trans.k]];
getP[g_, delta_] := Module[{deltatranspose = Transpose[delta, {1}]},
   deltatranspose.g.delta];

Now I need to use them for different values of the variable kappa, ranging from 1 to n. Even if the passages of the different calculations performed are a bit complicated to follow, the crucial point is that whenever kappa changes, not only getA, but also all the other functions, change. What I am really interested in is the values of the getW and getP functions. For every value of kappa from 1 to n, I should be able to store the values of getW and getP into two different "arrays" (actually, a matrix in the case of getW and a proper array in the case of getP). How can I do that, without using a For loop?

EDIT:

I am including the full code, where I specify what all the elements are, and I try to implement the suggestions received in the comments:

getA[kappa_] := 
  Table[2*Cos[(2*Pi/n)*(Abs[j - i])*kappa], {i, n}, {j, n}];
getF[csi_, a_, b_] := 
  Module[{csiInv = Inverse[csi]}, .5 Tr[csiInv.a.csiInv.b]];
getG[csi_, f_, a_] := Module[{csiInv = Inverse[csi]}, 
   csiInv.a.csiInv/2];
getE[g_, k_] := Module[{kinv = Inverse[k]},
   Transpose[kinv].g.kinv];
getW[k_, a_, e_] := Module[{ktrans = Transpose[k]},
   Tr[k.a.ktrans.e]];
getV[csi_, e_, e2_, k_] := 
  Module[{ktrans = Transpose[k], e2trans = Transpose[e2]},
   2*Tr[csi.ktrans.e.k.csi.ktrans.e2trans.k]];
getP[g_, delta_] := Module[{deltatranspose = Transpose[delta, {1}]},
   deltatranspose.g.delta];

n = L = 4;
sigma = 3;
nyquist = n/2 + 1;
sampling = 8;
mu = 0.0;

powerspectrum[i_] := 
  Piecewise[{{0, i == 0}, {Exp[-(2*Pi*i*sigma/L)^2], 
     0 < i <= n/2}, {Exp[-(2*Pi*(n - i)*sigma/L)^2], n/2 < i <= n}}];
pts = Table[powerspectrum[i], {i, 0, n}] ;

inverse = InverseFourier[pts];

func[inverse_] := 
  Module[{n = Length[inverse], tup}, 
   tup = Cases[
     Tuples[Range[n], 2], {i_, j_} /; Abs[i - j] < (n - 1)/2];
   SparseArray[
    Thread[tup -> (inverse[[Abs[#1 - #2] + 1]] & @@@ tup)], {n, n}]];
CSI = func[inverse];

Kmatrix = 
  Table[((3.0*70.0*70.0*0.3)/(2.0*300000.0*300000.0))*((j + 1)*(i + 
        2 - (j + 1)))*(1.0 + (70.0/300000.0)*(j + 1)), {i, 0, 
    n - 1}, {j, 0, n - 1}];

leftREAL = 
  Table[RandomVariate[
    NormalDistribution[mu, Exp[-(2*Pi*k*sigma/L)^2]]], {k, n/2}];
rightREAL = Reverse[leftREAL] /. {x_, y_} -> {n - x, y};
fullREAL = Join[{0.0}, Most[leftREAL], rightREAL];

leftIMAGINARY = 
  Table[RandomVariate[
    NormalDistribution[mu, Exp[-(2*Pi*k*sigma/L)^2]]], {k, n/2 - 1}];
rightIMAGINARY = -Reverse[leftIMAGINARY] /. {x_, y_} -> {n - x, y};
fullIMAGINARY = Join[{0.0}, leftIMAGINARY, {0.0}, rightIMAGINARY];

fullfield = fullREAL + I*fullIMAGINARY;

fieldconfiguration = InverseFourier[fullfield];

Wab = Table[
   getE[getG[CSI, getF[CSI, getA[alpha], getA[beta]], getA[alpha], 
     Kmatrix]], {alpha, n}, {beta, n}] // MatrixForm
Pa = Table[
   getP[getG[CSI, getA[alpha], getA[alpha]], 
    fieldconfiguration], {alpha, n}] // MatrixForm

I am now getting the following problem:

Dot::dotsh: Tensors {{1.39036 +0.225877 I,0.046518 -0.874121 I,-0.573392+0.0884825 I,0.046518 +0.30145 I,0.272321 -0.137394 I},<<3>>,{0.272321 -0.137394 I,0.046518 +0.30145 I,-0.573392+0.0884825 I,0.046518 -0.874121 I,1.39036 +0.225877 I}} and {{2,0,-2,0},{0,2,0,-2},{-2,0,2,0},{0,-2,0,2}} have incompatible shapes. 

which means that I am multiplying matrices or arrays that have not got the same dimension, but I cannot see where I am making this mistake. Does anyone see where I am going wrong?

$\endgroup$
  • $\begingroup$ How are they related? How do the other functions change when kappa changes? $\endgroup$ – Marius Ladegård Meyer Apr 13 '15 at 14:38
  • $\begingroup$ @MariusLadegårdMeyer every time that you see an 'a' matrix in the definition of the functions, that is a matrix obtained through the geta[kappa_definition]. So, this creates a chain: if I change 'kappa', then 'a' changes and so on... $\endgroup$ – johnhenry Apr 13 '15 at 14:44
  • 1
    $\begingroup$ Ok, let's take getW for instance. What is the first argument k? The second argument is getA[kappa], and the third is getE[getG[csi,getF[csi,getA[kappa],b],getA[kappa],k] What is csi, and what is b? If you have all these, and you want a table of all kappa-values from 1 to n, simply use Table. $\endgroup$ – Marius Ladegård Meyer Apr 13 '15 at 14:53
  • 1
    $\begingroup$ CSI is 5x5 and getA[1] is 4x4 so the dot product in getF fails $\endgroup$ – Simon Woods Apr 13 '15 at 17:54
  • 1
    $\begingroup$ Just trace it back through the code - CSI is 5x5 because inverse is length 5. inverse is length 5 because pts is length 5. Maybe the table should run from 0 to n-1 instead of 0 to n? I haven't tried to understand what you're doing so I don't know. $\endgroup$ – Simon Woods Apr 13 '15 at 18:05
1
$\begingroup$

If you know the kappa values ahead of time (eg. {1,2,3,4,5} ) and the dimension of your getA matrix is 3, for example

kappa = Range[5];
n=3;

then use a Pure Function and Map it to the List of kappa values

getA[#] & /@ kappa
getW[IdentityMatrix[3], getA[#], IdentityMatrix[3]] & /@ kappa

{
{{2, -1, -1}, {-1, 2, -1}, {-1, -1, 2}},
{{2, -1, -1}, {-1, 2, -1}, {-1, -1, 2}},
{{2, 2, 2}, {2, 2, 2}, {2, 2, 2}},
{{2, -1, -1}, {-1, 2, -1}, {-1, -1, 2}},
{{2, -1, -1}, {-1, 2, -1}, {-1, -1, 2}}
}

{6, 6, 6, 6, 6}

If you want this to happen dynamically, you'll have to use a different approach, but given your mention of a procedural For loop, I suspect that's not your intention.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.