Log In | Get Help   
Home My Page Projects Code Snippets Project Openings Mareframe
Summary Activity Forums Tracker Lists Tasks Docs Surveys News SCM Files
[mareframe] Diff of /trunk/gadget/hooke.cc
[mareframe] / trunk / gadget / hooke.cc Repository:
ViewVC logotype

Diff of /trunk/gadget/hooke.cc

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 11, Thu Jul 23 19:00:38 2015 UTC revision 12, Fri Jul 24 18:36:24 2015 UTC
# Line 188  Line 188 
188    
189  /* given a point, look for a better one nearby, one coord at a time */  /* given a point, look for a better one nearby, one coord at a time */
190  #ifndef NO_OPENMP  #ifndef NO_OPENMP
191    /*
192     * function bestBeraby parallelized with OpenMP
193     * · 2 threads per coord to parallelize the calculation of +delta/-delta
194     * · parallelize the calculation of the best nearby of the coord
195     */
196  double OptInfoHooke::bestNearbyOMP(DoubleVector& delta, DoubleVector& point, double prevbest, IntVector& param) {  double OptInfoHooke::bestNearbyOMP(DoubleVector& delta, DoubleVector& point, double prevbest, IntVector& param) {
 //cout << "beastNearbyINI" << endl;  
197    double minf;//, ftmp;    double minf;//, ftmp;
198    int i, j, k;    int i, j, k;
199    DoubleVector z(point);    DoubleVector z(point);
# Line 215  Line 219 
219    }    }
220    
221    while ( i < nvars) {    while ( i < nvars) {
 //        for (int k=0;k<point.Size(); k++)  
 //                cout << z[k] << " " ;  
 //          cout << endl;  
222            if ((i + paral_tokens -1) >= nvars)            if ((i + paral_tokens -1) >= nvars)
223                    paral_tokens = nvars - i;                    paral_tokens = nvars - i;
 //cout << "i:" << i << endl;  
224            omp_set_dynamic(0);            omp_set_dynamic(0);
225            omp_set_nested(1);            omp_set_nested(1);
226  #pragma omp parallel for num_threads(paral_tokens) private(k)  #pragma omp parallel for num_threads(paral_tokens) private(k)
227            for (j = 0; j < paral_tokens; ++j) {            for (j = 0; j < paral_tokens; ++j) {
 //                cout << "thr_for:" << omp_get_num_threads()<<endl;  
 //                cout << "J:" << j << endl;  
228                    storage[j].z = z;                    storage[j].z = z;
229                    storage[j].delta = delta;                    storage[j].delta = delta;
230                    DoubleVector v1(z);                    DoubleVector v1(z);
# Line 237  Line 235 
235    
236  #pragma omp parallel sections num_threads(2)  #pragma omp parallel sections num_threads(2)
237                    {                    {
 //                        cout << "thr_sec:" << omp_get_num_threads()<<endl;  
238    #pragma omp section    #pragma omp section
239                            {                            {
240                    storage[j].ftmp = EcoSystems[j]->SimulateAndUpdate(v1);                    storage[j].ftmp = EcoSystems[j]->SimulateAndUpdate(v1);
 //                cout << "->" << j << endl;  
241                            }                            }
242    #pragma omp section    #pragma omp section
243                            {                            {
 //                                cout << "->" << j+paral_tokens << endl;  
244                    storage[j+paral_tokens].ftmp = EcoSystems[j+paral_tokens]->SimulateAndUpdate(v2);                    storage[j+paral_tokens].ftmp = EcoSystems[j+paral_tokens]->SimulateAndUpdate(v2);
245                            }                            }
246                    }                    }
 //cout << "-----" << endl;  
 //                cout <<i+j<< "-z["<< param[i+j]<<"]:" << (storage[j].z)[param[i+j]] << " - " << storage[j].ftmp << endl;  
247    
248                    if (storage[j].ftmp < minf) {                    if (storage[j].ftmp < minf) {
249                            storage[j].iters = 1;                            storage[j].iters = 1;
250                            storage[j].z[k] = v1[k];                            storage[j].z[k] = v1[k];
 //                        storage[j].delta[k] += delta[k];  
 //                      minf = ftmp;  
251                    } else {                    } else {
252                            storage[j].iters = 2;                            storage[j].iters = 2;
253                            storage[j].delta[k] = 0.0 - delta[k];                            storage[j].delta[k] = 0.0 - delta[k];
254                            if (storage[j+paral_tokens].ftmp < minf) {                            if (storage[j+paral_tokens].ftmp < minf) {
 //                                minf = ftmp;  
255                                    storage[j].ftmp = storage[j+paral_tokens].ftmp;                                    storage[j].ftmp = storage[j+paral_tokens].ftmp;
256                                    storage[j].z[k] = v2[k];                                    storage[j].z[k] = v2[k];
257                            }                            }
 //                        else  
 //                                storage[j].z[k] = point[k];  
258                    }                    }
259            }            }
260    
 //        cout << i <<"-z["<< param[i]<<"]:" <<z[param[i]] << " - " ;  
261            for (j = 0; j < paral_tokens; ++j) {            for (j = 0; j < paral_tokens; ++j) {
262                    i++;                    i++;
263                    iters += storage[j].iters;                    iters += storage[j].iters;
# Line 278  Line 265 
265                            minf = storage[j].ftmp;                            minf = storage[j].ftmp;
266                            z = storage[j].z;                            z = storage[j].z;
267                            delta = storage[j].delta;                            delta = storage[j].delta;
 //                        cout << storage[j].ftmp;  
268                            break;                            break;
269                    }                    }
270            }            }
271    }    }
272    
 //  omp_set_num_threads(numThr);  
   
273    for (i = 0; i < nvars; ++i)    for (i = 0; i < nvars; ++i)
274      point[i] = z[i];      point[i] = z[i];
275    return minf;    return minf;
# Line 337  Line 321 
321      bestx[i] = x[i];      bestx[i] = x[i];
322      trialx[i] = x[i];      trialx[i] = x[i];
323      param[i] = i;      param[i] = i;
     //FIXME rand_r  
324      delta[i] = ((2 * (rand() % 2)) - 1) * rho;  //JMB - randomise the sign      delta[i] = ((2 * (rand() % 2)) - 1) * rho;  //JMB - randomise the sign
325    }    }
326    
# Line 371  Line 354 
354      /* randomize the order of the parameters once in a while */      /* randomize the order of the parameters once in a while */
355      rchange = 0;      rchange = 0;
356      while (rchange < nvars) {      while (rchange < nvars) {
         //FIXME rand_r  
357        rnumber = rand() % nvars;        rnumber = rand() % nvars;
358        rcheck = 1;        rcheck = 1;
359        for (i = 0; i < rchange; i++)        for (i = 0; i < rchange; i++)
# Line 401  Line 383 
383  #ifdef NO_OPENMP  #ifdef NO_OPENMP
384      iters = EcoSystem->getFuncEval() - offset;      iters = EcoSystem->getFuncEval() - offset;
385  #endif  #endif
 //    cout << "newf:" << newf << "-" << iters<<endl;  
386      if (iters > hookeiter) {      if (iters > hookeiter) {
387        handle.logMessage(LOGINFO, "\nStopping Hooke & Jeeves optimisation algorithm\n");        handle.logMessage(LOGINFO, "\nStopping Hooke & Jeeves optimisation algorithm\n");
388        handle.logMessage(LOGINFO, "The optimisation stopped after", iters, "function evaluations");        handle.logMessage(LOGINFO, "The optimisation stopped after", iters, "function evaluations");
# Line 499  Line 480 
480  #ifdef NO_OPENMP  #ifdef NO_OPENMP
481        iters = EcoSystem->getFuncEval() - offset;        iters = EcoSystem->getFuncEval() - offset;
482  #endif  #endif
 //      cout << "newf:" << newf << "-" << iters<<endl;  
483        if (iters > hookeiter) {        if (iters > hookeiter) {
484          handle.logMessage(LOGINFO, "\nStopping Hooke & Jeeves optimisation algorithm\n");          handle.logMessage(LOGINFO, "\nStopping Hooke & Jeeves optimisation algorithm\n");
485          handle.logMessage(LOGINFO, "The optimisation stopped after", iters, "function evaluations");          handle.logMessage(LOGINFO, "The optimisation stopped after", iters, "function evaluations");
# Line 552  Line 532 
532    }    }
533  }  }
534    
535    /* Functions to perform the parallelization of the algorithm of HJ with OpenMP*/
536  #ifdef GADGET_OPENMP  #ifdef GADGET_OPENMP
   //TODO doc  
537  double OptInfoHooke::bestNearbyOMP2(DoubleVector& delta, DoubleVector& point, double prevbest, IntVector& param) {  double OptInfoHooke::bestNearbyOMP2(DoubleVector& delta, DoubleVector& point, double prevbest, IntVector& param) {
538  //cout << "beastNearbyINI" << endl;    double minf;
   double minf;//, ftmp;  
539    int i, j, k, ii;    int i, j, k, ii;
540    DoubleVector z(point);    DoubleVector z(point);
541    int bestId = 0;    int bestId = 0;
# Line 569  Line 548 
548    
549    minf = prevbest;    minf = prevbest;
550    
   
551    int paral_tokens, numThr, nvars = point.Size();    int paral_tokens, numThr, nvars = point.Size();
552    numThr = omp_get_max_threads ( );    numThr = omp_get_max_threads ( );
553    
# Line 583  Line 561 
561    for (ii=0; ii< paral_tokens; ii++) {    for (ii=0; ii< paral_tokens; ii++) {
562            i = 0;            i = 0;
563            while ( i < nvars) {            while ( i < nvars) {
         //        for (int k=0;k<point.Size(); k++)  
         //                cout << z[k] << " " ;  
         //          cout << endl;  
564                    if ((i + paral_tokens -1) >= nvars)                    if ((i + paral_tokens -1) >= nvars)
565                            paral_tokens = nvars - i;                            paral_tokens = nvars - i;
         //cout << "i:" << i << endl;  
566                    omp_set_dynamic(0);                    omp_set_dynamic(0);
567                    omp_set_nested(1);                    omp_set_nested(1);
568          #pragma omp parallel for num_threads(paral_tokens) private(k)          #pragma omp parallel for num_threads(paral_tokens) private(k)
569                    for (j = 0; j < paral_tokens; ++j) {                    for (j = 0; j < paral_tokens; ++j) {
         //                cout << "thr_for:" << omp_get_num_threads()<<endl;  
         //                cout << "J:" << j << endl;  
570                            storage[j].z = z;                            storage[j].z = z;
571                            storage[j].delta = delta;                            storage[j].delta = delta;
572                            DoubleVector v1(z);                            DoubleVector v1(z);
# Line 605  Line 577 
577    
578          #pragma omp parallel sections num_threads(2)          #pragma omp parallel sections num_threads(2)
579                            {                            {
         //                        cout << "thr_sec:" << omp_get_num_threads()<<endl;  
580            #pragma omp section            #pragma omp section
581                                    {                                    {
582                            storage[j].ftmp = EcoSystems[j]->SimulateAndUpdate(v1);                            storage[j].ftmp = EcoSystems[j]->SimulateAndUpdate(v1);
         //                cout << "->" << j << endl;  
583                                    }                                    }
584            #pragma omp section            #pragma omp section
585                                    {                                    {
         //                                cout << "->" << j+paral_tokens << endl;  
586                            storage[j+paral_tokens].ftmp = EcoSystems[j+paral_tokens]->SimulateAndUpdate(v2);                            storage[j+paral_tokens].ftmp = EcoSystems[j+paral_tokens]->SimulateAndUpdate(v2);
587                                    }                                    }
588                            }                            }
         //cout << "-----" << endl;  
         //                cout <<i+j<< "-z["<< param[i+j]<<"]:" << (storage[j].z)[param[i+j]] << " - " << storage[j].ftmp << endl;  
   
589                            if (storage[j].ftmp < minf) {                            if (storage[j].ftmp < minf) {
590                                    storage[j].iters = 1;                                    storage[j].iters = 1;
591                                    storage[j].z[k] = v1[k];                                    storage[j].z[k] = v1[k];
         //                        storage[j].delta[k] += delta[k];  
         //                      minf = ftmp;  
592                            } else {                            } else {
593                                    storage[j].iters = 2;                                    storage[j].iters = 2;
594                                    storage[j].delta[k] = 0.0 - delta[k];                                    storage[j].delta[k] = 0.0 - delta[k];
595                                    if (storage[j+paral_tokens].ftmp < minf) {                                    if (storage[j+paral_tokens].ftmp < minf) {
         //                                minf = ftmp;  
596                                            storage[j].ftmp = storage[j+paral_tokens].ftmp;                                            storage[j].ftmp = storage[j+paral_tokens].ftmp;
597                                            storage[j].z[k] = v2[k];                                            storage[j].z[k] = v2[k];
598                                    }                                    }
599                                    else iters += 2;                                    else iters += 2;
         //                                storage[j].z[k] = point[k];  
600                            }                            }
601                    }                    }
602    
         //        cout << i <<"-z["<< param[i]<<"]:" <<z[param[i]] << " - " ;  
603                    bestId = 0;                    bestId = 0;
604                    for (j = 1; j < paral_tokens; ++j) {                    for (j = 1; j < paral_tokens; ++j) {
605                            if (storage[j].ftmp < storage[bestId].ftmp)                            if (storage[j].ftmp < storage[bestId].ftmp)
# Line 649  Line 610 
610                            minf = storage[bestId].ftmp;                            minf = storage[bestId].ftmp;
611                            z = storage[bestId].z;                            z = storage[bestId].z;
612                            delta = storage[bestId].delta;                            delta = storage[bestId].delta;
         //                        cout << storage[j].ftmp;  
   
613                    }                    }
614    
615                    i += paral_tokens;                    i += paral_tokens;
616            }            }
617          }          }
618    
 //  omp_set_num_threads(numThr);  
619    delete[] storage;    delete[] storage;
620    for (i = 0; i < nvars; ++i)    for (i = 0; i < nvars; ++i)
621      point[i] = z[i];      point[i] = z[i];
# Line 709  Line 667 
667               bestx[i] = x[i];               bestx[i] = x[i];
668               trialx[i] = x[i];               trialx[i] = x[i];
669               param[i] = i;               param[i] = i;
              //FIXME rand_r  
670               delta[i] = ((2 * (rand() % 2)) - 1) * rho;  //JMB - randomise the sign               delta[i] = ((2 * (rand() % 2)) - 1) * rho;  //JMB - randomise the sign
671             }             }
672    
# Line 743  Line 700 
700               /* randomize the order of the parameters once in a while */               /* randomize the order of the parameters once in a while */
701               rchange = 0;               rchange = 0;
702               while (rchange < nvars) {               while (rchange < nvars) {
                 //FIXME rand_r  
703                 rnumber = rand() % nvars;                 rnumber = rand() % nvars;
704                 rcheck = 1;                 rcheck = 1;
705                 for (i = 0; i < rchange; i++)                 for (i = 0; i < rchange; i++)
# Line 773  Line 729 
729           #ifdef NO_OPENMP           #ifdef NO_OPENMP
730               iters = EcoSystem->getFuncEval() - offset;               iters = EcoSystem->getFuncEval() - offset;
731           #endif           #endif
          //    cout << "newf:" << newf << "-" << iters<<endl;  
732               if (iters > hookeiter) {               if (iters > hookeiter) {
733                 handle.logMessage(LOGINFO, "\nStopping Hooke & Jeeves optimisation algorithm\n");                 handle.logMessage(LOGINFO, "\nStopping Hooke & Jeeves optimisation algorithm\n");
734                 handle.logMessage(LOGINFO, "The optimisation stopped after", iters, "function evaluations");                 handle.logMessage(LOGINFO, "The optimisation stopped after", iters, "function evaluations");
# Line 786  Line 741 
741                 for (i = 0; i < nvars; i++)                 for (i = 0; i < nvars; i++)
742                   bestx[i] = trialx[i] * init[i];                   bestx[i] = trialx[i] * init[i];
743                 EcoSystem->storeVariables(score, bestx);                 EcoSystem->storeVariables(score, bestx);
 //             EcoSystem->writeOptValuesOMP();  
744                 return;                 return;
745               }               }
746    
# Line 872  Line 826 
826           #ifdef NO_OPENMP           #ifdef NO_OPENMP
827                 iters = EcoSystem->getFuncEval() - offset;                 iters = EcoSystem->getFuncEval() - offset;
828           #endif           #endif
          //      cout << "newf:" << newf << "-" << iters<<endl;  
829                 if (iters > hookeiter) {                 if (iters > hookeiter) {
830                   handle.logMessage(LOGINFO, "\nStopping Hooke & Jeeves optimisation algorithm\n");                   handle.logMessage(LOGINFO, "\nStopping Hooke & Jeeves optimisation algorithm\n");
831                   handle.logMessage(LOGINFO, "The optimisation stopped after", iters, "function evaluations");                   handle.logMessage(LOGINFO, "The optimisation stopped after", iters, "function evaluations");
# Line 885  Line 838 
838                   for (i = 0; i < nvars; i++)                   for (i = 0; i < nvars; i++)
839                     bestx[i] = trialx[i] * init[i];                     bestx[i] = trialx[i] * init[i];
840                   EcoSystem->storeVariables(score, bestx);                   EcoSystem->storeVariables(score, bestx);
 //               EcoSystem->writeOptValuesOMP();  
841                   return;                   return;
842                 }                 }
843               }               }
# Line 916  Line 868 
868                 score = bestf;                 score = bestf;
869                 handle.logMessage(LOGINFO, "\nHooke & Jeeves finished with a likelihood score of", score);                 handle.logMessage(LOGINFO, "\nHooke & Jeeves finished with a likelihood score of", score);
870                 EcoSystem->storeVariables(bestf, bestx);                 EcoSystem->storeVariables(bestf, bestx);
 //             EcoSystem->writeOptValuesOMP();  
871                 return;                 return;
872               }               }
873    
Legend:
Removed from v.11  
changed lines
  Added in v.12
root@forge.cesga.es
 ViewVC Help
Powered by ViewVC 1.0.0  

Powered By FusionForge