tree_projection.cpp
Go to the documentation of this file.
00001 /*
00002  *    This file is part of ACADO Toolkit.
00003  *
00004  *    ACADO Toolkit -- A Toolkit for Automatic Control and Dynamic Optimization.
00005  *    Copyright (C) 2008-2014 by Boris Houska, Hans Joachim Ferreau,
00006  *    Milan Vukov, Rien Quirynen, KU Leuven.
00007  *    Developed within the Optimization in Engineering Center (OPTEC)
00008  *    under supervision of Moritz Diehl. All rights reserved.
00009  *
00010  *    ACADO Toolkit is free software; you can redistribute it and/or
00011  *    modify it under the terms of the GNU Lesser General Public
00012  *    License as published by the Free Software Foundation; either
00013  *    version 3 of the License, or (at your option) any later version.
00014  *
00015  *    ACADO Toolkit is distributed in the hope that it will be useful,
00016  *    but WITHOUT ANY WARRANTY; without even the implied warranty of
00017  *    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00018  *    Lesser General Public License for more details.
00019  *
00020  *    You should have received a copy of the GNU Lesser General Public
00021  *    License along with ACADO Toolkit; if not, write to the Free Software
00022  *    Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA
00023  *
00024  */
00025 
00026 
00027 
00035 #include <acado/utils/acado_utils.hpp>
00036 #include <acado/symbolic_operator/symbolic_operator.hpp>
00037 #include <acado/symbolic_expression/symbolic_expression.hpp>
00038 #include <acado/symbolic_expression/constraint_component.hpp>
00039 
00040 
00041 BEGIN_NAMESPACE_ACADO
00042 
00043 
00044 
00045 int TreeProjection::count = 0;
00046 
00047 TreeProjection::TreeProjection( )
00048                :Projection(){
00049 
00050     variableType   = VT_INTERMEDIATE_STATE ;
00051     vIndex         = 0                     ;
00052     variableIndex  = 0                     ;
00053     argument       = 0                     ;
00054     ne             = NE_ZERO               ;
00055 }
00056 
00057 
00058 TreeProjection::TreeProjection( const std::string &name_ )
00059                :Projection( name_ ){
00060 
00061     variableType   = VT_INTERMEDIATE_STATE ;
00062     vIndex         = 0                     ;
00063     variableIndex  = 0                     ;
00064     argument       = 0                     ;
00065     ne             = NE_ZERO               ;
00066 }
00067 
00068 
00069 
00070 TreeProjection::TreeProjection( const TreeProjection &arg )
00071                :Projection(){
00072 
00073     copy(arg);
00074 
00075     if( arg.argument == 0 ){
00076         argument = 0;
00077     }
00078     else{
00079         argument = arg.argument;
00080         argument->nCount++;
00081     }
00082 
00083     ne = arg.ne;
00084 }
00085 
00086 
00087 TreeProjection::~TreeProjection(){
00088  
00089     if( argument != 0 ){
00090 
00091         if( argument->nCount == 0 ){
00092             delete argument;
00093             argument = 0;
00094         }
00095         else{
00096             argument->nCount--;
00097         }
00098     }
00099 }
00100 
00101 
00102 
00103 Operator& TreeProjection::operator=( const Operator &arg ){
00104 
00105     if( this != &arg ){
00106 
00107         if( argument != 0 ){
00108             if( argument->nCount == 0 ){
00109                 delete argument;
00110                 argument = 0;
00111             }
00112             else{
00113                 argument->nCount--;
00114             }
00115         }
00116 
00117         Operator *arg_tmp = arg.clone();
00118         TreeProjection *tp = dynamic_cast<TreeProjection *>(arg_tmp);
00119         Projection *p = dynamic_cast<Projection *>(arg_tmp);
00120         if( tp != 0 ) {
00121                 // special case: argument is a treeprojection
00122                 if( tp->argument != 0 ) {
00123                         argument = tp->argument->clone();
00124                         ne = argument->isOneOrZero();
00125                 }
00126                 else {
00127                         argument = 0;
00128                         ne = NE_NEITHER_ONE_NOR_ZERO;
00129                 }
00130                         copy(*tp);
00131 
00132         }
00133         else {
00134                 if( p != 0 ) {
00135                         // special case: argument is a projection
00136                         argument = 0;
00137                         ne = NE_NEITHER_ONE_NOR_ZERO;
00138                         copy(*p);
00139                 }
00140                 else {
00141                         // no special case: create a new treeprojection
00142                         argument = arg.clone() ;
00143                         vIndex   = count++;
00144                         variableIndex  = vIndex ;
00145 
00146                         curvature      = CT_UNKNOWN; // argument->getCurvature();
00147                         monotonicity   = MT_UNKNOWN; // argument->getMonotonicity();
00148 
00149                         ne = argument->isOneOrZero();
00150 
00151                         if( curvature == CT_CONSTANT )
00152                                 scale = argument->getValue();
00153                 }
00154         }
00155         delete arg_tmp;
00156     }
00157 
00158     return *this;
00159 }
00160 
00161 
00162 Operator& TreeProjection::operator=( const Expression &arg ){
00163 
00164         ASSERT( arg.getDim() == 1 );
00165 
00166         if( argument != 0 ){
00167                 if( argument->nCount == 0 ){
00168                         delete argument;
00169                         argument = 0;
00170                 }
00171                 else{
00172                         argument->nCount--;
00173                 }
00174         }
00175 
00176         argument = arg.getOperatorClone(0);
00177 
00178         vIndex         = count++;
00179         variableIndex  = vIndex ;
00180 
00181         curvature      = CT_UNKNOWN; // argument->getCurvature();
00182         monotonicity   = MT_UNKNOWN; // argument->getMonotonicity();
00183 
00184         ne = argument->isOneOrZero();
00185 
00186         if( curvature == CT_CONSTANT )
00187                 scale = argument->getValue();
00188 
00189         return *this;
00190 }
00191 
00192 
00193 Operator& TreeProjection::operator=( const double& arg ){
00194 
00195     scale = arg;
00196 
00197     Expression tmp( arg );
00198     return this->operator=( tmp );
00199 }
00200 
00201 
00202 TreeProjection* TreeProjection::clone() const{
00203 
00204     return new TreeProjection( *this );
00205 }
00206 
00207 
00208 TreeProjection* TreeProjection::cloneTreeProjection() const{
00209 
00210     return new TreeProjection( *this );
00211 }
00212 
00213 
00214 
00215 Operator* TreeProjection::ADforwardProtected( int dim,
00216                                                    VariableType *varType,
00217                                                    int *component,
00218                                                    Operator **seed,
00219                                                    int &nNewIS,
00220                                                    TreeProjection ***newIS ){
00221 
00222         if (argument == 0) {
00223                 return Projection::ADforwardProtected( dim, varType, component, seed, nNewIS, newIS );
00224         }
00225         else {
00226 
00227     int run1 = 0;
00228 
00229     while( run1 < dim ){
00230 
00231         if( varType[run1] == variableType && component[run1] == vIndex ){
00232             return seed[run1]->clone();
00233         }
00234         run1++;
00235     }
00236 
00237     if( vIndex >= nNewIS ){
00238 
00239         *newIS = (TreeProjection**)realloc(*newIS,(vIndex+1)*sizeof(TreeProjection*));
00240 
00241         for( run1 = nNewIS; run1 < vIndex + 1; run1++ )
00242              newIS[0][run1] = 0;
00243 
00244         nNewIS = vIndex+1;
00245     }
00246 
00247     if( newIS[0][vIndex] != 0 ){
00248 
00249         return newIS[0][vIndex]->clone();
00250     }
00251 
00252     Operator *tmp = argument->AD_forward(dim,varType,component,seed,nNewIS,newIS);
00253 
00254     newIS[0][vIndex] = new TreeProjection();
00255 
00256     newIS[0][vIndex]->operator=(*tmp);
00257 
00258     newIS[0][vIndex]->setCurvature( CT_UNKNOWN );
00259 
00260     delete tmp;
00261     return newIS[0][vIndex]->clone();
00262         }
00263 }
00264 
00265 
00266 
00267 returnValue TreeProjection::ADbackwardProtected( int           dim      , 
00268                                         VariableType *varType  , 
00269                                         int          *component, 
00270                                         Operator     *seed     , 
00271                                         Operator    **df       , 
00272                                         int           &nNewIS  , 
00273                                         TreeProjection ***newIS   ){
00274 
00275         if (argument == 0) {
00276                 return Projection::ADbackwardProtected( dim, varType, component, seed, df, nNewIS, newIS );
00277         }
00278         else {
00279 
00280     int run1;
00281 
00282     if( (vIndex+1)*dim-1 >= nNewIS ){
00283 
00284         *newIS = (TreeProjection**)realloc(*newIS,((vIndex+1)*dim)*sizeof(TreeProjection*));
00285 
00286         for( run1 = nNewIS; run1 < (vIndex+1)*dim; run1++ )
00287              newIS[0][run1] = 0;
00288 
00289         nNewIS = (vIndex+1)*dim;
00290     }
00291     
00292     if( newIS[0][vIndex*dim] == 0 ){
00293       
00294         Operator **results = new Operator*[dim];
00295         for( run1 = 0; run1 < dim; run1++ ){
00296             results[run1] = new DoubleConstant(0.0,NE_ZERO);
00297         }
00298         Operator *aux = new DoubleConstant(1.0,NE_ONE);
00299         argument->AD_backward(dim,varType,component, aux, results, nNewIS, newIS );
00300         for( run1 = 0; run1 < dim; run1++ ){  
00301               newIS[0][vIndex*dim+run1] = new TreeProjection();
00302               newIS[0][vIndex*dim+run1]->operator=(*results[run1]);
00303               newIS[0][vIndex*dim+run1]->setCurvature( CT_UNKNOWN );
00304               delete results[run1];
00305         }
00306         delete[] results;
00307     }
00308 
00309     if( seed->isOneOrZero() != NE_ZERO ){
00310       
00311         for( run1 = 0; run1 < dim; run1++ ){
00312           
00313             Operator *tmp = df[run1]->clone();
00314             delete df[run1];
00315             
00316             if( seed->isOneOrZero() == NE_ONE ){
00317                 df[run1] = myAdd( newIS[0][vIndex*dim+run1], tmp );
00318             }
00319             else{
00320                 Operator *projTmp = myProd(newIS[0][vIndex*dim+run1],seed);
00321                 df[run1] = myAdd( projTmp, tmp );
00322                 delete projTmp;
00323             }
00324             delete tmp;
00325         }
00326     }
00327     delete seed;
00328 
00329     return SUCCESSFUL_RETURN;
00330         }
00331 }
00332 
00333 
00334 
00335 returnValue TreeProjection::ADsymmetricProtected( int            dim       , 
00336                                         VariableType  *varType   , 
00337                                         int           *component , 
00338                                         Operator      *l         , 
00339                                         Operator     **S         , 
00340                                         int            dimS      , 
00341                                         Operator     **dfS       , 
00342                                         Operator     **ldf       , 
00343                                         Operator     **H         , 
00344                                         int            &nNewLIS  , 
00345                                         TreeProjection ***newLIS , 
00346                                         int            &nNewSIS  , 
00347                                         TreeProjection ***newSIS , 
00348                                         int            &nNewHIS  , 
00349                                         TreeProjection ***newHIS    ){
00350   
00351   
00352         if (argument == 0) {
00353                 return Projection::ADsymmetricProtected( dim, varType, component, l, S, dimS, dfS, ldf, H, nNewLIS, newLIS, nNewSIS, newSIS, nNewHIS, newHIS );
00354         }
00355         else {
00356 
00357         int run1;
00358 
00359         if( (vIndex+1)*dim-1 >= nNewLIS ){
00360 
00361                 *newLIS = (TreeProjection**)realloc(*newLIS,((vIndex+1)*dim)*sizeof(TreeProjection*));
00362                 for( run1 = nNewLIS; run1 < (vIndex+1)*dim; run1++ )
00363                         newLIS[0][run1] = 0;
00364                 nNewLIS = (vIndex+1)*dim;
00365         }
00366 
00367         if( (vIndex+1)*dimS-1 >= nNewSIS ){
00368 
00369                 *newSIS = (TreeProjection**)realloc(*newSIS,((vIndex+1)*dimS)*sizeof(TreeProjection*));
00370                 for( run1 = nNewSIS; run1 < (vIndex+1)*dimS; run1++ )
00371                         newSIS[0][run1] = 0;
00372                 nNewSIS = (vIndex+1)*dimS;
00373         }
00374 
00375         if( (vIndex+1)*dimS*dimS-1 >= nNewHIS ){
00376 
00377                 *newHIS = (TreeProjection**)realloc(*newHIS,((vIndex+1)*dimS*dimS)*sizeof(TreeProjection*));
00378                 for( run1 = nNewHIS; run1 < (vIndex+1)*dimS*dimS; run1++ )
00379                         newHIS[0][run1] = 0;
00380                 nNewHIS = (vIndex+1)*dimS*dimS;
00381         }
00382   
00383   // ============================================================================
00384 
00385         if( newLIS[0][vIndex*dim] == 0 ){
00386 
00387                 Operator **lres = new Operator*[dim];
00388                 for( run1 = 0; run1 < dim; run1++ ){
00389                         lres[run1] = new DoubleConstant(0.0,NE_ZERO);
00390                 }
00391                 Operator *aux = new DoubleConstant(1.0,NE_ONE);
00392 
00393                 Operator **Sres = new Operator*[dimS];
00394                 for( run1 = 0; run1 < dimS; run1++ ){
00395                         Sres[run1] = new DoubleConstant(0.0,NE_ZERO);
00396                 }
00397 
00398                 Operator **Hres = new Operator*[dimS*dimS];
00399                 for( run1 = 0; run1 < dimS*dimS; run1++ ){
00400                         Hres[run1] = new DoubleConstant(0.0,NE_ZERO);
00401                 }
00402 
00403                 argument->AD_symmetric( dim, varType, component, aux, S, dimS,
00404                                 Sres, lres, Hres,
00405                                 nNewLIS, newLIS, nNewSIS, newSIS, nNewHIS, newHIS );
00406 
00407 
00408                 for( run1 = 0; run1 < dim; run1++ ){
00409                         newLIS[0][vIndex*dim+run1] = new TreeProjection();
00410                         newLIS[0][vIndex*dim+run1]->operator=(*lres[run1]);
00411                         newLIS[0][vIndex*dim+run1]->setCurvature( CT_UNKNOWN );
00412                         delete lres[run1];
00413                 }
00414                 delete[] lres;
00415 
00416                 for( run1 = 0; run1 < dimS; run1++ ){
00417                         newSIS[0][vIndex*dimS+run1] = new TreeProjection();
00418                         newSIS[0][vIndex*dimS+run1]->operator=(*Sres[run1]);
00419                         newSIS[0][vIndex*dimS+run1]->setCurvature( CT_UNKNOWN );
00420                         delete Sres[run1];
00421                 }
00422                 delete[] Sres;
00423 
00424                 for( run1 = 0; run1 < dimS*dimS; run1++ ){
00425                         newHIS[0][vIndex*dimS*dimS+run1] = new TreeProjection();
00426                         newHIS[0][vIndex*dimS*dimS+run1]->operator=(*Hres[run1]);
00427                         newHIS[0][vIndex*dimS*dimS+run1]->setCurvature( CT_UNKNOWN );
00428                         delete Hres[run1];
00429                 }
00430                 delete[] Hres;
00431         }
00432 
00433   // ============================================================================
00434   
00435         if( l->isOneOrZero() != NE_ZERO ){
00436 
00437                 for( run1 = 0; run1 < dim; run1++ ){
00438 
00439                         Operator *tmp = ldf[run1]->clone();
00440                         delete ldf[run1];
00441 
00442                         if( l->isOneOrZero() == NE_ONE ){
00443                                 ldf[run1] = myAdd( newLIS[0][vIndex*dim+run1], tmp );
00444                         }
00445                         else{
00446                                 Operator *projTmp = myProd( newLIS[0][vIndex*dim+run1], l );
00447                                 ldf[run1] = myAdd( projTmp, tmp );
00448                                 delete projTmp;
00449                         }
00450                         delete tmp;
00451                 }
00452         }
00453 
00454         for( run1 = 0; run1 < dimS; run1++ ){
00455                 delete dfS[run1];
00456                 dfS[run1] = newSIS[0][vIndex*dimS+run1]->clone();
00457         }
00458 
00459         for( run1 = 0; run1 < dimS*dimS; run1++ ){
00460 
00461                 delete H[run1];
00462 
00463                 if( l->isOneOrZero() == NE_ONE ){
00464                         H[run1] = newHIS[0][vIndex*dimS*dimS+run1]->clone();
00465                 }
00466                 else{
00467                         H[run1] = myProd(newHIS[0][vIndex*dimS*dimS+run1],l);
00468                 }
00469         }
00470 
00471         delete l;
00472         return SUCCESSFUL_RETURN;
00473         }
00474 }
00475   
00476 
00477   
00478 returnValue TreeProjection::loadIndices( SymbolicIndexList *indexList ){
00479 
00480     returnValue returnvalue = SUCCESSFUL_RETURN;
00481 
00482     if( argument == 0 ){
00483         return Projection::loadIndices( indexList );
00484     }
00485 
00486     if( indexList->addNewElement( VT_INTERMEDIATE_STATE, vIndex ) == BT_TRUE ){
00487 
00488         returnvalue = argument->loadIndices( indexList );
00489 
00490         indexList->addOperatorPointer( argument, vIndex );
00491     }
00492 
00493     if (name.empty())
00494     {
00495         std::stringstream ss;
00496         ss << "a" << "[" << vIndex << "]";
00497         name = ss.str();
00498     }
00499 
00500     return returnvalue;
00501 }
00502 
00503 
00504 BooleanType TreeProjection::isDependingOn( int dim,
00505                                                 VariableType *varType,
00506                                                 int *component,
00507                                                 BooleanType   *implicit_dep ){
00508 
00509         if( argument == 0 ){
00510                 return Projection::isDependingOn( dim, varType, component, implicit_dep );
00511         }
00512         return implicit_dep[vIndex];
00513 }
00514 
00515 
00516 BooleanType TreeProjection::isLinearIn( int dim,
00517                                              VariableType *varType,
00518                                              int *component,
00519                                              BooleanType   *implicit_dep ){
00520 
00521         if( argument == 0 ){
00522                 return Projection::isLinearIn( dim, varType, component, implicit_dep );
00523         }
00524     return implicit_dep[vIndex];
00525 }
00526 
00527 
00528 BooleanType TreeProjection::isPolynomialIn( int dim,
00529                                                  VariableType *varType,
00530                                                  int *component,
00531                                                  BooleanType   *implicit_dep ){
00532 
00533         if( argument == 0 ){
00534                 return Projection::isPolynomialIn( dim, varType, component, implicit_dep );
00535         }
00536     return implicit_dep[vIndex];
00537 }
00538 
00539 
00540 BooleanType TreeProjection::isRationalIn( int dim,
00541                                                VariableType *varType,
00542                                                int *component,
00543                                                BooleanType   *implicit_dep ){
00544 
00545         if( argument == 0 ){
00546                 return Projection::isRationalIn( dim, varType, component, implicit_dep );
00547         }
00548         return implicit_dep[vIndex];
00549 }
00550 
00551 
00552 returnValue TreeProjection::clearStaticCounters(){
00553 
00554     count     = 0;
00555     return SUCCESSFUL_RETURN;
00556 }
00557 
00558 
00559 Operator* TreeProjection::getArgument() const{
00560 
00561     if( argument != 0 ) return argument->clone();
00562     return 0;
00563 }
00564 
00565 
00566 NeutralElement TreeProjection::isOneOrZero() const{
00567 
00568     return ne;
00569 }
00570 
00571 
00572 
00573 void TreeProjection::copy( const Projection &arg ){
00574 
00575     Projection::copy(arg);
00576 }
00577 
00578 
00579 Operator* TreeProjection::passArgument() const{
00580 
00581     return argument;
00582 }
00583 
00584 
00585 returnValue TreeProjection::setVariableExportName(      const VariableType &_type,
00586                                                                                                         const std::vector< std::string >& _name
00587                                                                                                         )
00588 {
00589         if (argument != 0 && argument->getName() == ON_POWER_INT)
00590                 argument->setVariableExportName(_type, _name);
00591 
00592         return Projection::setVariableExportName(_type, _name);
00593 }
00594 
00595 
00596 BooleanType TreeProjection::isTrivial() const {
00597 
00598         return argument != 0;
00599 }
00600 
00601 
00602 returnValue TreeProjection::initDerivative() {
00603 
00604         if( !initialized && argument != 0 ) {
00605                 initialized = BT_TRUE;
00606                 return argument->initDerivative();
00607         }
00608         else {
00609                 return SUCCESSFUL_RETURN;
00610         }
00611 }
00612 
00613 
00614 CLOSE_NAMESPACE_ACADO
00615 
00616 // end of file.


acado
Author(s): Milan Vukov, Rien Quirynen
autogenerated on Sat Jun 8 2019 19:40:00