I just copy-pasted your script into my code. I changed only the name of the finite element space and I wrote “us” instead of “v” (as I already have “v” defined as radial velocity) and M which is the equivalent of “DistA” matrix.
fespace XXXXFES(Th,[P2, P2, P2, P1]);
fespace XFES (Th, P2);XXXXFES< complex >[int] [ EigenVecu,EigenVecv,EigenVecw,EigenVecp ](nev); //array to store eigenvector
XXXXFES< complex >[int] [ AdjEigenVecu,AdjEigenVecv,AdjEigenVecw,AdjEigenVecp ](nev); //array to store adjoint eigenvector
int i, k, Adjk;
k = EPSSolve(M, J, values=EigenVal, vectors=EigenVecu, sparams = …); //Solve eigenvalue problem /direct
Adjk = EPSSolve(AdjM, AdjJ, values=AdjEigenVal, vectors=AdjEigenVecu, sparams = …); //Solve eigenvalue problem /adjoint
real[int] D;
for (i=0; i<k; i++)
{complex wrong = int2d(Th)(EigenVecu[i]*AdjEigenVecu[i]+EigenVecv[i]*AdjEigenVecv[i]);
cout << mpirank << " " << wrong << endl;
XFES< complex > scaledVecu = EigenVecu[i];
for[i, us : M.D] scaledVecu[i] = us;
XFES< complex > scaledVecv = EigenVecv[i][];
for[i, us : M.D] scaledVecv[][i] = us;
wrong = int2d(Th)(scaledVecuAdjEigenVecu[i]+scaledVecvAdjEigenVecv[i]);
complex reduce;
mpiAllReduce(wrong, reduce, mpiCommWorld, mpiSUM);
if(mpirank == 0)
cout << reduce << endl;}
As soon as it reach *for[i, us : M.D] scaledVecu[i] = us; I get the error from above. I guess the way I define the array to store eigenvector causes the problem with “array type”.