kaskr / adcomp

AD computation with Template Model Builder (TMB)
Other
178 stars 80 forks source link

Creating an Array from a Matrix #183

Closed Njoselson closed 8 years ago

Njoselson commented 8 years ago

Hi,

I am trying to set all the elements of a matrix equal to an equally sized array.

  matrix<Type> PPMmyc(105,16);     //(1,105,0,aplus);   //l,a
  PPMmyc.setZero();
  array<Type> PPMmyca(105,16);
  PPMmyca = PPMmyc.array(); 

This compiles correctly, but gives a fatal error on running the objective function.

I see from the Eigen documentation that there is an ArrayWrapper class. http://eigen.tuxfamily.org/dox/classEigen_1_1ArrayWrapper.html I am unclear if this is the same functionality as m1.array() but regardless, I can't get the syntax to compile correctly.

Also I see the Map Class seems to do something similar as well. https://eigen.tuxfamily.org/dox/group__TutorialMapClass.html

Is there a TMB specific function for changing from a Matrix to an Array?

Nathaniel

Njoselson commented 8 years ago

In the array definition dox, it says that setting an array equal to a matrix is allowed. Line 149 https://github.com/kaskr/adcomp/blob/7ff8bb3a1d048a5d52c9fb1271605bbf0a37cbe6/TMB/inst/include/tmbutils/array.hpp However, my session crashes with a fatal error. When I run the debugger the output is as follows:

TMB has received an error from Eigen. The following condition was not met:
rows() == other.rows() && cols() == other.cols()
Please check your matrix-vector bounds etc., or run your program through a debugger.
Process 2912 stopped
* thread #1: tid = 0x27f49, 0x00007fff9ad30f06 libsystem_kernel.dylib`__pthread_kill + 10, queue = 'com.apple.main-thread', stop reason = signal SIGABRT
    frame #0: 0x00007fff9ad30f06 libsystem_kernel.dylib`__pthread_kill + 10
libsystem_kernel.dylib`__pthread_kill:
->  0x7fff9ad30f06 <+10>: jae    0x7fff9ad30f10            ; <+20>
    0x7fff9ad30f08 <+12>: movq   %rax, %rdi
    0x7fff9ad30f0b <+15>: jmp    0x7fff9ad2b7cd            ; cerror_nocancel
    0x7fff9ad30f10 <+20>: retq 
(lldb) bt
* thread #1: tid = 0x27f49, 0x00007fff9ad30f06 libsystem_kernel.dylib`__pthread_kill + 10, queue = 'com.apple.main-thread', stop reason = signal SIGABRT
  * frame #0: 0x00007fff9ad30f06 libsystem_kernel.dylib`__pthread_kill + 10
    frame #1: 0x00007fffa098a4ec libsystem_pthread.dylib`pthread_kill + 90
    frame #2: 0x00007fff9bd316e7 libsystem_c.dylib`abort + 129
    frame #3: 0x000000010a812e9f Test2.so`Eigen::Map<Eigen::Array<double, -1, 1, 0, -1, 1>, 0, Eigen::Stride<0, 0> >& Eigen::DenseBase<Eigen::Map<Eigen::Array<double, -1, 1, 0, -1, 1>, 0, Eigen::Stride<0, 0> > >::lazyAssign<Eigen::Matrix<double, -1, -1, 0, -1, -1> >(this=<unavailable>, other=<unavailable>) + 1311 at Assign.h:505
    frame #4: 0x000000010a803058 Test2.so`objective_function<double>::operator()() [inlined] Eigen::internal::assign_selector<Eigen::Map<Eigen::Array<double, -1, 1, 0, -1, 1>, 0, Eigen::Stride<0, 0> >, Eigen::Matrix<double, -1, -1, 0, -1, -1>, false, false>::run(Eigen::Map<Eigen::Array<double, -1, 1, 0, -1, 1>, 0, Eigen::Stride<0, 0> >&, Eigen::Matrix<double, -1, -1, 0, -1, -1> const&) + 2296 at Assign.h:527
    frame #5: 0x000000010a803045 Test2.so`objective_function<double>::operator()() [inlined] Eigen::Map<Eigen::Array<double, -1, 1, 0, -1, 1>, 0, Eigen::Stride<0, 0> >& Eigen::DenseBase<Eigen::Map<Eigen::Array<double, -1, 1, 0, -1, 1>, 0, Eigen::Stride<0, 0> > >::operator=<Eigen::Matrix<double, -1, -1, 0, -1, -1> >(this=<unavailable>, other=<unavailable>) at Assign.h:552
    frame #6: 0x000000010a803045 Test2.so`objective_function<double>::operator()() [inlined] Eigen::Map<Eigen::Array<double, -1, 1, 0, -1, 1>, 0, Eigen::Stride<0, 0> >& Eigen::Map<Eigen::Array<double, -1, 1, 0, -1, 1>, 0, Eigen::Stride<0, 0> >::operator=<Eigen::Matrix<double, -1, -1, 0, -1, -1> >(other=<unavailable>) at Map.h:170
    frame #7: 0x000000010a803045 Test2.so`objective_function<double>::operator()() [inlined] tmbutils::array<double> tmbutils::array<double>::operator=<tmbutils::matrix<double> >(this=<unavailable>) at array.hpp:270
    frame #8: 0x000000010a803045 Test2.so`objective_function<double>::operator(this=<unavailable>)() + 2277 at Test2.cpp:25
    frame #9: 0x000000010a828fe7 Test2.so`::getParameterOrder(data=<unavailable>, parameters=<unavailable>, report=<unavailable>) + 103 at tmb_core.hpp:1150

It references line 25, which is:

22  matrix<Type> PPMmyc(105,16);     //(1,105,0,aplus);   //l,a
23  PPMmyc.setZero();
24  array<Type> PPMmyca(105,16);
25  PPMmyca = PPMmyc;

What am I missing which is causing the error? All my bounds seem to match.

Sorry if I am asking obvious questions.

Nathaniel

skaug commented 8 years ago

Regarding your first message that "but gives a fatal error on running the objective function". This does not happen to me. I attached a modified example that shows that your code produce the correct result. issue_183.zip

However, I admit that array/matrix is a bit confusion, and we should work on clarifying it. Some comments have been added to the "documentation" branch, and these will find their way into master in the future.

If you after working on this see ways to improve the docs, please make a pull request to branch documentation.

kaskr commented 8 years ago

@Njoselson Note some array changes took place in v1.6.6 - are you using a recent version?: https://github.com/kaskr/adcomp/blob/master/TMB/NEWS

Njoselson commented 8 years ago

Hi, Thank you for your response.

My version of TMB is the 1.7.0. I just re-cloned the TMB folder and did a new make install from the command line.

Running your example I still get a fatal error. My debug report is attached:

TMB has received an error from Eigen. The following condition was not met:
rows() == other.rows() && cols() == other.cols()
Please check your matrix-vector bounds etc., or run your program through a debugger.
Process 1463 stopped
* thread #1: tid = 0xf50a, 0x00007fff9315df06 libsystem_kernel.dylib`__pthread_kill + 10, queue = 'com.apple.main-thread', stop reason = signal SIGABRT
    frame #0: 0x00007fff9315df06 libsystem_kernel.dylib`__pthread_kill + 10
libsystem_kernel.dylib`__pthread_kill:
->  0x7fff9315df06 <+10>: jae    0x7fff9315df10            ; <+20>
    0x7fff9315df08 <+12>: movq   %rax, %rdi
    0x7fff9315df0b <+15>: jmp    0x7fff931587cd            ; cerror_nocancel
    0x7fff9315df10 <+20>: retq   
(lldb) bt
* thread #1: tid = 0xf50a, 0x00007fff9315df06 libsystem_kernel.dylib`__pthread_kill + 10, queue = 'com.apple.main-thread', stop reason = signal SIGABRT
  * frame #0: 0x00007fff9315df06 libsystem_kernel.dylib`__pthread_kill + 10
    frame #1: 0x00007fff98db74ec libsystem_pthread.dylib`pthread_kill + 90
    frame #2: 0x00007fff9415e6e7 libsystem_c.dylib`abort + 129
    frame #3: 0x000000010b158bc0 arrays_183.so`Eigen::Map<Eigen::Array<double, -1, 1, 0, -1, 1>, 0, Eigen::Stride<0, 0> >& Eigen::DenseBase<Eigen::Map<Eigen::Array<double, -1, 1, 0, -1, 1>, 0, Eigen::Stride<0, 0> > >::lazyAssign<Eigen::ArrayWrapper<Eigen::Matrix<double, -1, -1, 0, -1, -1> > >(this=<unavailable>, other=<unavailable>) + 832 at Assign.h:505
    frame #4: 0x000000010b14e380 arrays_183.so`objective_function<double>::operator()() [inlined] Eigen::internal::assign_selector<Eigen::Map<Eigen::Array<double, -1, 1, 0, -1, 1>, 0, Eigen::Stride<0, 0> >, Eigen::ArrayWrapper<Eigen::Matrix<double, -1, -1, 0, -1, -1> >, false, false>::run(Eigen::Map<Eigen::Array<double, -1, 1, 0, -1, 1>, 0, Eigen::Stride<0, 0> >&, Eigen::ArrayWrapper<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&) + 384 at Assign.h:527
    frame #5: 0x000000010b14e370 arrays_183.so`objective_function<double>::operator()() [inlined] Eigen::Map<Eigen::Array<double, -1, 1, 0, -1, 1>, 0, Eigen::Stride<0, 0> >& Eigen::DenseBase<Eigen::Map<Eigen::Array<double, -1, 1, 0, -1, 1>, 0, Eigen::Stride<0, 0> > >::operator=<Eigen::ArrayWrapper<Eigen::Matrix<double, -1, -1, 0, -1, -1> > >(this=<unavailable>, other=<unavailable>) at Assign.h:552
    frame #6: 0x000000010b14e370 arrays_183.so`objective_function<double>::operator()() [inlined] Eigen::Map<Eigen::Array<double, -1, 1, 0, -1, 1>, 0, Eigen::Stride<0, 0> >& Eigen::Map<Eigen::Array<double, -1, 1, 0, -1, 1>, 0, Eigen::Stride<0, 0> >::operator=<Eigen::ArrayWrapper<Eigen::Matrix<double, -1, -1, 0, -1, -1> > >(other=<unavailable>) at Map.h:170
    frame #7: 0x000000010b14e370 arrays_183.so`objective_function<double>::operator()() [inlined] tmbutils::array<double> tmbutils::array<double>::operator=<Eigen::ArrayWrapper<Eigen::Matrix<double, -1, -1, 0, -1, -1> > >(this=<unavailable>) at array.hpp:270
    frame #8: 0x000000010b14e370 arrays_183.so`objective_function<double>::operator(this=<unavailable>)() + 368 at arrays_183.cpp:18
    frame #9: 0x000000010b171e07 arrays_183.so`getParameterOrder + 103

It references line 18

PPMmyca = PPMmyc.array(); 

Perhaps it is a problem with my eigen, or my clang? I am running Mac OS X 10.11.4

Thanks!

Nathaniel

Njoselson commented 8 years ago

Is there a way to re-install as version 1.6.6 to see if that helps?

So neither of you get a fatal error by running this example?

kaskr commented 8 years ago

Can't reproduce this problem with v1.7.0. Have tried both gcc and clang.

My guess: you somehow didn't update TMB correctly (the error from the debugger indicates an old version). Could it be caused by multiple installations of R? TMB? Or perhaps an old precompiled TMB?

Njoselson commented 8 years ago

I will delete all previous TMB Folders and delete R version 3.1, and then tell you what happens.

Njoselson commented 8 years ago

I deleted, re-cloned and reinstalled TMB but the problem persists. I am using MRO version 3.2.3, which I just installed a few weeks ago.

Is there any way to determine if the update to version 1.7.0 failed?

Njoselson commented 8 years ago

Reinstalled the software on my computer and reinstalled everything and now the compilation works!

Thanks so much for all your help!