             vec[ieq][0] = mat[0][0]*pdsoln[0] + mat[0][1]*pdsoln[1];
             vec[ieq][1] = mat[1][0]*pdsoln[0] + mat[1][1]*pdsoln[1];
 #endif
-            pdsoln += NDIM;
+            pdsoln += (NDIM+1)*NDIM;
+        };
+        pten = (double(*)[NDIM])(pidsoln+NDIM);
+        for (it=0; it<NDIM; it++) {
+            for (jt=0; jt<NDIM; jt++) {
+                vmt[it][jt] = mat[it][0]*pten[0][jt] + mat[it][1]*pten[1][jt]
+#if NDIM == 3
+                            + mat[it][2]*pten[2][jt]
+#endif
+                ;
+            };
+        };
+        for (it=0; it<NDIM; it++) {
+            for (jt=0; jt<NDIM; jt++) {
+                vec[it+1][jt] = vmt[it][0]*mvt[0][jt] + vmt[it][1]*mvt[1][jt]
+#if NDIM == 3
+                              + vmt[it][2]*mvt[2][jt]
+#endif
+                ;
+            };
         };
         // set wall condition in the rotated coordinate;
         vec[0][0] = -vec[0][0];

    double len;
#endif
    // arrays.
    double mat[NDIM][NDIM], matinv[NDIM][NDIM];
    double mom[NDIM];
    // iterators.
    int ifc, icl, jcl;

    int *pfacn, *pfccls;
    double *pidsoln, *pjdsoln, *pdsoln;
    double *pfcnml;
    double (*pten)[NDIM];
#if NDIM == 3
    int *pfcnds;
    double *pndcrd, *pfccnd;
