investment.cpp

Go to the documentation of this file.
00001 
00016 #include <iomanip>
00017 #include "flopc.hpp"
00018 
00019 // Change these two lines to use a different solver
00020 #include "OsiClpSolverInterface.hpp"
00021 #define OSI_SOLVER_INTERFACE OsiClpSolverInterface
00022 
00023 using namespace std;
00024 using namespace flopc;
00025 
00026 
00028 
00033 class StageNode {
00034         public:
00035                 MP_set ASSETS;                
00036                 MP_index a;                   
00037                 MP_variable x;                
00038                 MP_expression objFunction;    
00039 
00040                 StageNode *ptParent;          
00041                 vector<StageNode *> children; 
00042                 int nodeNmb;                  
00043                 double prob;                  
00044 
00045                 StageNode(StageNode *ptPred, const int nmbAssets, const int nodeN,
00046                           const double condProb)
00047                 : ASSETS(nmbAssets), x(ASSETS), ptParent(ptPred), nodeNmb(nodeN),
00048                   prob(condProb)
00049                 {
00050                         if (ptParent != NULL) {
00051                                 prob *= ptParent->prob;             // Compute the total probability
00052                                 ptParent->children.push_back(this); // Register with the parent
00053                         }
00054                 }
00055                 virtual ~StageNode(){}
00056 
00057         protected:
00059 
00061                 virtual void make_obj_function_() {
00062                         assert (children.size() > 0 && "Leaves should never call this function");
00063                         // Create the objective recursively for all descendants of the node
00064                         for (int i=0; i < (int) children.size(); i++) {
00065                                 children[i]->make_obj_function_();
00066                         }
00067                         // No objective at non-leaf nodes, so the objective value at this node
00068                         // and all its descendants is a sum of obj. values of the children.
00069                         objFunction = children[0]->objFunction;
00070                         for (int i=1; i < (int) children.size(); i++) {
00071                                 objFunction = objFunction + children[i]->objFunction;
00072                         }
00073                 }
00074 };
00075 
00077 class RootNode : public StageNode {
00078         public:
00079                 MP_constraint initialBudget; 
00080 
00081                 RootNode(const int nmbAssets, const double initWealth)
00082                 : StageNode(NULL, nmbAssets, 0, 1.0)
00083                 {
00084                         initialBudget() = sum(ASSETS, x(ASSETS)) == initWealth;
00085                 }
00086 
00088 
00090                 void make_objective_function() {
00091                         make_obj_function_();
00092                 }
00093 };
00094 
00096 class MidStageNode : public StageNode {
00097         public:
00098                 MP_constraint cashFlowBalance; 
00099                 MP_data Return;                
00100 
00102 
00108                 MidStageNode(StageNode *ptPred, const int nodeN, const double condProb,
00109                              double *ptRetVect)
00110                 : StageNode(ptPred, ptPred->ASSETS.size(), nodeN, condProb),
00111                   Return(ptRetVect, ASSETS)
00112                 {
00113                         // This shows the use of MP_index in a formula
00114                         cashFlowBalance = sum(ASSETS(a), ptParent->x(a) * Return(a))
00115                                           == sum(ASSETS(a), x(a));
00116                 }
00117                 virtual ~MidStageNode(){}
00118 };
00119 
00121 class LeafNode : public StageNode {
00122         public:
00123                 MP_variable w;              
00124                 MP_variable y;              
00125                 MP_constraint finalBalance; 
00126                 MP_data Return;             
00127 
00129 
00133                 LeafNode(StageNode *ptPred, const int nodeN, const double condProb,
00134                          const double *ptRetVect, const double capTarget)
00135                 : StageNode(ptPred, ptPred->ASSETS.size(), nodeN, condProb), Return(ASSETS)
00136                 {
00137                         Return.value(ptRetVect); // Copy values from retVect to Return
00138                         // This shows a formula without using an additional MP_index
00139                         finalBalance = sum(ASSETS, ptParent->x(ASSETS) * Return(ASSETS))
00140                                        + w() - y() == capTarget;
00141                 }
00142 
00143         protected:
00145                 void make_obj_function_() {
00146                         objFunction = prob * (1.3 * w() - 1.1 * y());
00147                 }
00148 };
00149 
00150 
00151 // -------------------------------------------------------------------------- //
00152 
00154 class ScenTreeStruct {
00155         public:
00156                 int nmbNodes;   
00157                 int firstLeaf;  
00158 
00160                 ScenTreeStruct(const int nNodes, const int firstL)
00161                 : nmbNodes(nNodes), firstLeaf(firstL)
00162                 {}
00163                 virtual ~ScenTreeStruct(){}
00164 
00166 
00169                 virtual int get_parent(int n) const = 0;
00170 };
00171 
00173 class BinTreeStruct : public ScenTreeStruct {
00174         public:
00176                 BinTreeStruct(const int nmbStages)
00177                 : ScenTreeStruct(pow(2.0, nmbStages) - 1, pow(2.0, nmbStages-1) - 1) {}
00178 
00179                 int get_parent(int n) const {
00180                         return (n-1) / 2;    // This gives: get_parent(0) = 0
00181                 }
00182 };
00183 
00184 // -------------------------------------------------------------------------- //
00185 
00186 int main()
00187 {
00188         // DATA - This would be normally read from some external file
00189         // binary scenario tree with 4 stages: 15 nodes, firstLeaf = 7
00190         BinTreeStruct scTree(4);
00191 
00192         // model parameters
00193         enum {stocks, bonds, nmbAssets}; // assets to invest to; sets nmbAssets = 2
00194         double InitBudget = 55;          // initial budget
00195         double CapTarget = 80;           // capital target
00196 
00197         // vector of returns at the 14 non-root nodes
00198         double retData[][nmbAssets] = {
00199                 {1.25, 1.14},
00200                 {1.06, 1.16},
00201                 {1.21, 1.17},
00202                 {1.07, 1.12},
00203                 {1.15, 1.18},
00204                 {1.06, 1.12},
00205                 {1.26, 1.13},
00206                 {1.07, 1.14},
00207                 {1.25, 1.15},
00208                 {1.06, 1.12},
00209                 {1.05, 1.17},
00210                 {1.06, 1.15},
00211                 {1.05, 1.14},
00212                 {1.06, 1.12}
00213         };
00214 
00215         // ------------------------------------------------------------------------ //
00216         // Initialize the object for the model
00217         MP_model &model = MP_model::getDefaultModel();
00218         model.setSolver(new OSI_SOLVER_INTERFACE);
00219         #ifdef NDEBUG
00220                 model.silent();  // non-debug, i.e. release version -> less output
00221         #else
00222                 model.verbose(); // debug version -> show more output
00223         #endif
00224 
00225         // Create the scenario tree
00226         vector<StageNode *> treeNodes;
00227         treeNodes.push_back(new RootNode(nmbAssets, InitBudget));
00228         for (int i = 1; i < scTree.nmbNodes; i++) {
00229                 int pred = scTree.get_parent(i);
00230                 if (i < scTree.firstLeaf) {
00231                         // Adding a mid-stage node
00232                         treeNodes.push_back(
00233                                 // Remember that this uses a shallow copy of retData,
00234                                 // i.e. that it links the constraints to the retData array!
00235                                 new MidStageNode(treeNodes[pred], i, 1.0 / 2, retData[i-1])
00236                         );
00237                 } else {
00238                         // Adding a leaf
00239                         treeNodes.push_back(
00240                                 // This constructor uses a deep copy, so it does not create a link.
00241                                 new LeafNode(treeNodes[pred], i, 1.0 / 2, retData[i-1], CapTarget)
00242                         );
00243                 }
00244         }
00245 
00246         #ifndef NDEBUG
00247         // Write some info about the scenario tree
00248         for (int n = 0; n < scTree.nmbNodes; n++) {
00249                 cout << "Node " << setw(2) << treeNodes[n]->nodeNmb << " has "
00250                      << "prob = " << setw(5) << treeNodes[n]->prob << " and "
00251                      << treeNodes[n]->children.size() << " children: ";
00252                 for (int cn = 0; cn < (int) treeNodes[n]->children.size(); cn++)
00253                         cout << setw(2) << treeNodes[n]->children[cn]->nodeNmb << " ";
00254                 cout << endl;
00255         }
00256         #endif
00257 
00258         // create a "shortcut object" for the root
00259         assert (dynamic_cast<RootNode *>(treeNodes[0]) != NULL
00260                                         && "The first node in treeNodes should be the Root");
00261         RootNode *ptRoot = dynamic_cast<RootNode *>(treeNodes[0]);
00262 
00263         ptRoot->make_objective_function();       // Create the objective func.
00264         model.setObjective(ptRoot->objFunction); // Set the objective
00265         model.attach();                          // Attach the model
00266         model->writeMps("investment");           // Write an MPS file
00267 
00268         // Solve the model
00269         model.minimize();
00270 
00271         // Report some results
00272         cout << endl << "Optimal objective value = " << model->getObjValue() << endl
00273                          << "First-stage x:";
00274         ptRoot->x.display();
00275 
00276         for (int n = 0; n < (int) treeNodes.size(); n++)
00277                 delete treeNodes[n]; // cleaning
00278 
00279         return 0;
00280 }

Generated on Mon Jan 26 10:52:55 2009 for investment by  doxygen 1.5.4