NDArray oct file gives error with 4 dimension

classic Classic list List threaded Threaded
3 messages Options
Reply | Threaded
Open this post in threaded view
|

NDArray oct file gives error with 4 dimension

niconeuman
Greetings,
I'm writing an oct-file for performing an operation between a 2D NDArray and
a 4D NDArray. Apparently trying to do

retMat(mu,nu) += D(ka,la)*gmnkl(mu,nu,ka,la);

does not compile. On the other hand, doing

retMat(mu,nu) += D(ka,ka)*gmnkl(mu,nu,ka);

compiles with no problem. Can anyone help me?
Best regards,
Nicolas

Below is the code for contractgD.cc,


#include <octave/oct.h>

DEFUN_DLD (contractgD, args, nargout,
           "This function takes as arguments D and gmnkl, and calculates the
matrix"
           "Jtemp which results from contracting gmnkl with D along
dimensions 3 and 4")
{

  NDArray D = args(0).array_value ();
  NDArray gmnkl = args(1).array_value ();

  const int lengthmu = gmnkl.dims()(0); //have to use it like this, not
NDArray.dims(i);
  const int lengthnu = gmnkl.dims()(1);
  const int lengthka = gmnkl.dims()(2);
  const int lengthla = gmnkl.dims()(3);

  octave_stdout << lengthmu << "\n";
  octave_stdout << lengthnu << "\n";
  octave_stdout << lengthka << "\n";
  octave_stdout << lengthla << "\n";
  octave_stdout << "gmnkl.ndims() = "  << "\n";
  octave_stdout << gmnkl.ndims() << "\n";
  octave_stdout << "gmnkl(1,1,1) = "  << "\n";
  octave_stdout << gmnkl(1,1,1) << "\n";
  octave_stdout << "gmnkl(1,1,1,1) = "  << "\n";
  octave_stdout << gmnkl(1,1,1,1) << "\n";
  dim_vector dvmn (lengthmu,lengthnu);

  NDArray retMat (dvmn);

  if (gmnkl.ndims() == 4){
  for(octave_idx_type mu = 0; mu < lengthmu; mu++){
    for(octave_idx_type nu = 0; nu < lengthnu; nu++){
      for(octave_idx_type ka = 0; ka < lengthka; ka++){
        for(octave_idx_type la = 0; la < lengthla; la++){
          retMat(mu,nu) += D(ka,la)*gmnkl(mu,nu,ka,la);
        }
      }
    }
  }
}

  octave_value_list retval (nargout);
  retval(0) = octave_value(retMat);
  return retval;
}



--
Sent from: https://octave.1599824.n4.nabble.com/Octave-General-f1599825.html


Reply | Threaded
Open this post in threaded view
|

Re: NDArray oct file gives error with 4 dimension

siko1056
On 2/7/20 6:52 AM, niconeuman wrote:

> Greetings,
> I'm writing an oct-file for performing an operation between a 2D NDArray and
> a 4D NDArray. Apparently trying to do
>
> retMat(mu,nu) += D(ka,la)*gmnkl(mu,nu,ka,la);
>
> does not compile. On the other hand, doing
>
> retMat(mu,nu) += D(ka,ka)*gmnkl(mu,nu,ka);
>
> compiles with no problem. Can anyone help me?
> Best regards,
> Nicolas
>
> Below is the code for contractgD.cc,
>
>
> #include <octave/oct.h>
>
> DEFUN_DLD (contractgD, args, nargout,
>            "This function takes as arguments D and gmnkl, and calculates the
> matrix"
>            "Jtemp which results from contracting gmnkl with D along
> dimensions 3 and 4")
> {
>
>   NDArray D = args(0).array_value ();
>   NDArray gmnkl = args(1).array_value ();
>
>   const int lengthmu = gmnkl.dims()(0); //have to use it like this, not
> NDArray.dims(i);
>   const int lengthnu = gmnkl.dims()(1);
>   const int lengthka = gmnkl.dims()(2);
>   const int lengthla = gmnkl.dims()(3);
>
>   octave_stdout << lengthmu << "\n";
>   octave_stdout << lengthnu << "\n";
>   octave_stdout << lengthka << "\n";
>   octave_stdout << lengthla << "\n";
>   octave_stdout << "gmnkl.ndims() = "  << "\n";
>   octave_stdout << gmnkl.ndims() << "\n";
>   octave_stdout << "gmnkl(1,1,1) = "  << "\n";
>   octave_stdout << gmnkl(1,1,1) << "\n";
>   octave_stdout << "gmnkl(1,1,1,1) = "  << "\n";
>   octave_stdout << gmnkl(1,1,1,1) << "\n";
>   dim_vector dvmn (lengthmu,lengthnu);
>
>   NDArray retMat (dvmn);
>
>   if (gmnkl.ndims() == 4){
>   for(octave_idx_type mu = 0; mu < lengthmu; mu++){
>     for(octave_idx_type nu = 0; nu < lengthnu; nu++){
>       for(octave_idx_type ka = 0; ka < lengthka; ka++){
>         for(octave_idx_type la = 0; la < lengthla; la++){
>           retMat(mu,nu) += D(ka,la)*gmnkl(mu,nu,ka,la);
>         }
>       }
>     }
>   }
> }
>
>   octave_value_list retval (nargout);
>   retval(0) = octave_value(retMat);
>   return retval;
> }
>
Unfortunately, the ease of ND-array handling of the Octave language does
not translate 1:1 to it's C++ API.  You get the "comfort" of ND-indices
by calling

   gmnkl(dim_vector (1,1,1,1).as_array ())

instead of

   gmnkl(1,1,1,1)

Same for

   retMat(mu,nu) += D(ka,la)*
      gmnkl(dim_vector (mu,nu,ka,la).as_array ());

Full program attached.  Tested in Octave 5.2.0 with

   mkoctfile contractgD.cc
   D = eye (2);
   gmnkl = eye (2, 2, 2, 2)
   contractgD (D, gmnkl)

As further suggestion, you should check for dimension at the very
beginning of your program to avoid bad memory access crashes, caused by
improper input arrays.

Read more at


https://octave.org/doxygen/5.2/d1/db1/classdim__vector.html#a0189767745d163089c5145b36adf24d4

https://octave.org/doxygen/5.2/d3/dd5/classNDArray.html#a803cb3dcd25b55bcc1cd5333f6373cb1

HTH,
Kai



contractgD.cc (1K) Download Attachment
Reply | Threaded
Open this post in threaded view
|

Re: NDArray oct file gives error with 4 dimension

niconeuman
Thank you very much Kai,
I solved it by transforming the 4D indices to a 1D index.
I will definitely add many checks at the beginning, because the code needs
to handle 1D to 4D cases.
Thanks again,
Best,
Nicolas




--
Sent from: https://octave.1599824.n4.nabble.com/Octave-General-f1599825.html