ADMB Documentation  11.1.2274
 All Classes Files Functions Variables Typedefs Friends Defines
fvar_m14.cpp
Go to the documentation of this file.
00001 /*
00002  * $Id: fvar_m14.cpp 1919 2014-04-22 22:02:01Z johnoel $
00003  *
00004  * Author: David Fournier
00005  * Copyright (c) 2008-2012 Regents of the University of California
00006  */
00011 // constructors, destructors and misc functions involving class prevariable
00012 
00013 #include "fvar.hpp"
00014 
00015 #ifdef __TURBOC__
00016   #pragma hdrstop
00017   #include <iostream.h>
00018 #endif
00019 
00020 #ifdef __ZTC__
00021   #include <iostream.hpp>
00022 #endif
00023 
00024 /*
00025  dvar_matrix  operator*(const dvar_matrix& m1, const dvar_matrix& m2 )
00026  {
00027    if (m1.colmin() != m2.rowmin() || m1.colmax() != m2.rowmax())
00028    {
00029      cerr << " Incompatible array bounds in dmatrix  operator * (const dmatrix& x, const dmatrix& m)\n";
00030      ad_exit(21);
00031    }
00032    //dmatrix cm1=value(m1);
00033    //dmatrix cm2=value(m2);
00034    dmatrix tmp(m1.rowmin(),m1.rowmax(), m2.colmin(), m2.colmax());
00035    double sum;
00036    double ** temp_col=(double **) malloc(m2.rowsize()*sizeof(double*));
00037    temp_col-=m2.rowmin();
00038 
00039 
00040    for (int j=m2.colmin(); j<=m2.colmax(); j++)
00041    {
00042 
00043      for (int k=m2.rowmin(); k<=m2.rowmax(); k++)
00044      {
00045        temp_col[k] = (double*) &m2.elem_value(k,j);
00046      }
00047 
00048      for (int i=m1.rowmin(); i<=m1.rowmax(); i++)
00049      {
00050        sum=0.0;
00051        const dvar_vector& temp_row = m1(i);
00052        for (int k=m1.colmin(); k<=m1.colmax(); k++)
00053        {
00054           sum+=temp_row.elem_value(k) * (*temp_col[k]);
00055          // sum+=temp_row(k) * cm2(k,j);
00056        }
00057        tmp(i,j)=sum;
00058      }
00059    }
00060 
00061 
00062    temp_col+=m2.rowmin();
00063    free ((char*)temp_col);
00064    dvar_matrix vtmp=nograd_assign(tmp);
00065    save_identifier_string("TEST1");
00066    m1.save_dvar_matrix_value();
00067    m1.save_dvar_matrix_position();
00068    m2.save_dvar_matrix_value();
00069    m2.save_dvar_matrix_position();
00070    vtmp.save_dvar_matrix_position();
00071    save_identifier_string("TEST6");
00072    gradient_structure::GRAD_STACK1->
00073             set_gradient_stack(dmdm_prod);
00074    return vtmp;
00075  }
00076 */
00077 
00082  dvar_matrix operator*(const dvar_matrix& m1, const dvar_matrix& m2)
00083  {
00084    if (m1.colmin() != m2.rowmin() || m1.colmax() != m2.rowmax())
00085    {
00086      cerr << " Incompatible array bounds in "
00087      "dmatrix operator*(const dmatrix& x, const dmatrix& m)\n";
00088      ad_exit(21);
00089    }
00090    //dmatrix cm1=value(m1);
00091    //dmatrix cm2=value(m2);
00092    dmatrix tmp(m1.rowmin(),m1.rowmax(), m2.colmin(), m2.colmax());
00093    double sum;
00094 
00095 
00096    for (int j=m2.colmin(); j<=m2.colmax(); j++)
00097    {
00098      dvector m2col=column_value(m2,j);
00099 
00100      for (int i=m1.rowmin(); i<=m1.rowmax(); i++)
00101      {
00102        sum=value(m1(i))*m2col;
00103        tmp(i,j)=sum;
00104      }
00105    }
00106 
00107    dvar_matrix vtmp=nograd_assign(tmp);
00108    save_identifier_string("TEST1");
00109    m1.save_dvar_matrix_value();
00110    m1.save_dvar_matrix_position();
00111    m2.save_dvar_matrix_value();
00112    m2.save_dvar_matrix_position();
00113    vtmp.save_dvar_matrix_position();
00114    save_identifier_string("TEST6");
00115    gradient_structure::GRAD_STACK1->
00116             set_gradient_stack(dmdm_prod);
00117    return vtmp;
00118  }
00119 
00124 void dmdm_prod(void)
00125 {
00126   verify_identifier_string("TEST6");
00127   dvar_matrix_position vpos=restore_dvar_matrix_position();
00128   dmatrix dftmp=restore_dvar_matrix_derivatives(vpos);
00129   dvar_matrix_position m2pos=restore_dvar_matrix_position();
00130   dmatrix cm2=restore_dvar_matrix_value(m2pos);
00131   dvar_matrix_position m1pos=restore_dvar_matrix_position();
00132   dmatrix cm1=restore_dvar_matrix_value(m1pos);
00133   verify_identifier_string("TEST1");
00134   dmatrix dfm1(m1pos);
00135   dmatrix dfm2(m2pos);
00136   double dfsum;
00137   dfm1.initialize();
00138   dfm2.initialize();
00139   for (int j=cm2.colmin(); j<=cm2.colmax(); j++)
00140   {
00141     for (int i=cm1.rowmin(); i<=cm1.rowmax(); i++)
00142     {
00143       //tmp.elem(i,j)=sum;
00144       dfsum=dftmp.elem(i,j);
00145       for (int k=cm1.colmin(); k<=cm1.colmax(); k++)
00146       {
00147         //sum+=cm1(i,k) * cm2(k,j);
00148         dfm1.elem(i,k)+=dfsum * cm2.elem(k,j);
00149         dfm2.elem(k,j)+=dfsum * cm1.elem(i,k);
00150       }
00151     }
00152   }
00153   dfm1.save_dmatrix_derivatives(m1pos);
00154   dfm2.save_dmatrix_derivatives(m2pos);
00155   // cout << "leaving dmdm_prod"<<endl;
00156 }