Min/Max

how to find max/min of this function

Page 1 of 1

11 Replies - 1913 Views - Last Post: 27 October 2010 - 03:19 AM Rate Topic: -----

#1 araven7  Icon User is offline

  • New D.I.C Head

Reputation: 0
  • View blog
  • Posts: 4
  • Joined: 26-October 10

Min/Max

Posted 26 October 2010 - 04:52 PM

One thing i need to do in order to start answering a given problem for this program is finding the maximum or minimum of: energy(y)
Can anyone help please? I would love a quick reply, I need this asap.
Thanks.

I know how to scan an array and test a variable eg: if min > array[i] set min to array[i], but i do not know how to get energy(y) into an array like array[i].

here's the program:

#include <iostream>  
#include <cmath>   
#include <algorithm> 
using namespace std;  
            
void rhs (double rho[], double f[]);   
void step_RK2(double h, double yin[], double yout[]);   
double energy (double y[]);    
          
main() {  
     
   double rhox,vx,tau;    
   double rhoy,vy;    
   double h,hnew;   
   double E0;         
     double y[4],yh[4],y2h[4];   
    int i;    
          
double safety=0.9,max_increase=2,max_decrease=5;    
double ratio,eps,eps0,atol,rtol;  
int accept,p,two_to_p;   
h=0.01;
atol=rtol=1e-3;  
tau=0;   
rhox=0;rhoy=1;   
vy=0.4;     
vx=0.4;        
      
/* the following lines initialise array y;  
         they also serve as a translation table  
        between vector elements and physical variables */    

y[0]=rhox; 
y[1]=rhoy;  
y[2]=vx;     
y[3]=vy;         
E0=energy(y);
     
       cout << tau << "\t" << rhox << "    \t"  << rhoy << "    \t"  << vx        
     << "    \t" << vy             << "    \t" << E0    
       << "    \t" << 0            << "    \t" << h       
     << "\n" ;         

while (tau < 300)  {          
        /* two smaller steps */    
     step_RK2(h,y ,yh);     
     step_RK2(h,yh,yh);        
        /* one double step */   
     step_RK2(2*h,y ,y2h);       
        /* error estimate for 2nd-order method*/     
     p=2;two_to_p=4;       
        ratio=0;   
        for(i=0;i<4;i++) {      
         eps=abs(y2h[i] - yh[i])/(two_to_p - 1);      
         eps0=atol+abs(yh[i])*rtol;    
   	 ratio=max(ratio,eps/eps0);      

                         }             

     /* step size correction */  
 
       hnew=h/exp( log(ratio)/(p+1) );         
            /* accept/reject the step */     
      if(ratio<1) {        
      if(accept==1) h=min(max_increase*h,safety*hnew);      

            for(i=0;i<4;i++)
              y[i]=(two_to_p*yh[i] - y2h[i])/(two_to_p - 1);     
              accept=1;       
              tau += 2*h;
          
                 } 

      else {          
       accept=0;           
       h=max(h/max_decrease,safety*hnew);    
           };   
       
       cout << tau << "\t" << y[0] << "    \t"  << y[1] << "    \t"         
         << y[2] << "    \t" << y[3]                 << "    \t" << energy(y)     
         << "    \t" << energy(y)-E0              << "    \t" << h           
      << "\n" ;             
        } ;      
  } ; 
     
void rhs (double y[], double force[]) {       
double rho_abs,rho3; 
     int i;          
rho_abs=sqrt(y[0]*y[0]+y[1]*y[1]);   
rho3=rho_abs*rho_abs*rho_abs;       
force[0] = y[2]; /* vx */   
force[1] = y[3]; /* vx */     
force[2] = -y[0]/rho3; /* -rhox/rho3 */   
force[3] = -y[1]/rho3; /* -rhoy/rho3 */      
       };            
void step_RK2(double h, double yin[], double yout[]) {     
double ytemp[4],k1[4],k2[4];     
int i;        
/* evaluate k1 */        
rhs(yin,k1);            
/* intermediate step */   
for(i=0;i<4;i++) ytemp[i] = yin[i] + h*k1[i];     
     /* evaluate k2 */    
       rhs(ytemp,k2);             
     /* update y */       
    for(i=0;i<4;i++) yout[i] = yin[i] + h*(k1[i]+k2[i])/2;   
          };         
double energy (double y[]) {         
return( (y[2]*y[2]+y[3]*y[3])/2 - 1/sqrt(y[0]*y[0]+y[1]*y[1]) );      
                           }

This post has been edited by Jayman: 26 October 2010 - 06:47 PM


Is This A Good Question/Topic? 0
  • +

Replies To: Min/Max

#2 Ollie9  Icon User is offline

  • D.I.C Head

Reputation: 6
  • View blog
  • Posts: 91
  • Joined: 10-May 10

Re: Min/Max

Posted 26 October 2010 - 04:53 PM

Put the code inside the [code] tags
Was This Post Helpful? 0
  • +
  • -

#3 araven7  Icon User is offline

  • New D.I.C Head

Reputation: 0
  • View blog
  • Posts: 4
  • Joined: 26-October 10

Re: Min/Max

Posted 26 October 2010 - 05:01 PM

View PostOllie9, on 26 October 2010 - 03:53 PM, said:

Put the code inside the [code] tags


fixed, I think.
Was This Post Helpful? 0
  • +
  • -

#4 janotte  Icon User is offline

  • code > sword
  • member icon

Reputation: 990
  • View blog
  • Posts: 5,141
  • Joined: 28-September 06

Re: Min/Max

Posted 26 October 2010 - 05:03 PM

Welcome to DIC!

As well as putting your code in code tags, like this :code: as Ollie9 has said also answer these questions.
( a ) Does your code compile?
( b ) Any errors or warnings? If there are then share them with us.
( c ) Is the program producing any output?
( d ) How is the actual output different to what you want / expect? Give details and, ideally, examples.
( e ) What have you already tried to fix it?
Was This Post Helpful? 0
  • +
  • -

#5 lfaivor21  Icon User is offline

  • New D.I.C Head

Reputation: 10
  • View blog
  • Posts: 29
  • Joined: 17-October 10

Re: Min/Max

Posted 26 October 2010 - 05:12 PM

Wouldn't think you'd need to put it in an array. Just test each energy value that comes out of energy inside your loop and set it to a min max variable if that value is bigger or smaller than the current values. If you don't mind me asking... what system is represented by this energy function? I see a .5*Velocity^2-1/r. Body under gravitational acceleration?

This post has been edited by lfaivor21: 26 October 2010 - 05:12 PM

Was This Post Helpful? 1
  • +
  • -

#6 araven7  Icon User is offline

  • New D.I.C Head

Reputation: 0
  • View blog
  • Posts: 4
  • Joined: 26-October 10

Re: Min/Max

Posted 26 October 2010 - 05:23 PM

View Postjanotte, on 26 October 2010 - 04:03 PM, said:

Welcome to DIC!

As well as putting your code in code tags, like this :code: as Ollie9 has said also answer these questions.
( a ) Does your code compile?
( b ) Any errors or warnings? If there are then share them with us.
( c ) Is the program producing any output?
( d ) How is the actual output different to what you want / expect? Give details and, ideally, examples.
( e ) What have you already tried to fix it?


a) yep it compiles as it is, i just want to add a few lines to it to find the max of energy(y)

b,) no errors/warnings.

c) yes it outputs several columns of data, each column corresponds to a set calculated variable. eg column 6 is all the values of energy(y)

d) i plan on just calculating the max/min value of energy(y) and have the code print this value at the very end of the output.
eg

tau - rhox - rhoy - vx - vy - energy(y) - energy(y)-E0 - h
0 0 1 0.4 0.4 0 0 0.023

columns of values

row of final values
Max of energy(y) = ??

e) i've tried adding in a loop to get:
for(i=0;i<n:i++){
array[i] = energy(y)}

and then I could just use a 'if min > array[i] set min to array[i]' type thing to find it, but it doesnt work, im sure its just i dont understand how to access the elements in energy(y).

View Postlfaivor21, on 26 October 2010 - 04:12 PM, said:

Wouldn't think you'd need to put it in an array. Just test each energy value that comes out of energy inside your loop and set it to a min max variable if that value is bigger or smaller than the current values. If you don't mind me asking... what system is represented by this energy function? I see a .5*Velocity^2-1/r. Body under gravitational acceleration?


yep, a kepler problem.

This post has been edited by araven7: 26 October 2010 - 05:21 PM

Was This Post Helpful? 0
  • +
  • -

#7 lfaivor21  Icon User is offline

  • New D.I.C Head

Reputation: 10
  • View blog
  • Posts: 29
  • Joined: 17-October 10

Re: Min/Max

Posted 26 October 2010 - 05:40 PM

Since gravity is a conservative vector field shouldn't the total energy be constant throughout? So why do you need max min except to check for errors? and double check the sign on your scalar potential because it seems you defined the force to be F = -1/r^3 (xi + yj) so the negative gradient of that potential should give you your force. But the one you have posted will give F = + 1/r^3(xi+yj) Might screw up the constant energy thing.

This post has been edited by lfaivor21: 26 October 2010 - 05:43 PM

Was This Post Helpful? 1
  • +
  • -

#8 janotte  Icon User is offline

  • code > sword
  • member icon

Reputation: 990
  • View blog
  • Posts: 5,141
  • Joined: 28-September 06

Re: Min/Max

Posted 26 October 2010 - 05:45 PM

View Postaraven7, on 27 October 2010 - 09:23 AM, said:

i dont understand how to access the elements in energy(y).


I'm not sure I understand what you are trying to say there.

energy() is a function.
You pass it an array and the values held in various elements of the array are used to calculate a return value.
So do you mean you are not sure you are extracting the values from the array to use in the calculation correctly?
I really don't understand what you are asking there.

At the moment you only call energy() once.
Don't you need to call it for all the sets of values you have? Something like that for loop you said you tried for loading an array.
Sorry I am just confused with what you are trying to do and why what you tried did not work.

BTW - Here's some feedback from the compiler
DIC.cpp:10: warning: ISO C++ forbids declaration of ‘main’ with no type
DIC.cpp: In function ‘void rhs(double*, double*)’:
DIC.cpp:88: warning: unused variable ‘i’


Was This Post Helpful? 0
  • +
  • -

#9 lfaivor21  Icon User is offline

  • New D.I.C Head

Reputation: 10
  • View blog
  • Posts: 29
  • Joined: 17-October 10

Re: Min/Max

Posted 26 October 2010 - 05:55 PM

he calls energy again here inside the while(tau<300) loop

		   cout << tau << "\t" << y[0] << "    \t"  << y[1] << "    \t"         
			<< y[2] << "    \t" << y[3]
			<< "    \t" << energy(y)     
			<< "    \t" << energy(y)-E0 
			<< "    \t" << h           
		        << "\n" ; 



Two ways you can do this... if you don't need to keep all the values returned by energy() as you vary y you can just test each new value to see if it's bigger or smaller than some max min and if it is then assign that new max min to max min. If you want to keep all the values assign each value to an element in an array and then loop through the array again with the same testing idea i laid out previously.
Was This Post Helpful? 1
  • +
  • -

#10 araven7  Icon User is offline

  • New D.I.C Head

Reputation: 0
  • View blog
  • Posts: 4
  • Joined: 26-October 10

Re: Min/Max

Posted 26 October 2010 - 06:03 PM

View Postlfaivor21, on 26 October 2010 - 04:55 PM, said:

he calls energy again here inside the while(tau<300) loop

		   cout << tau << "\t" << y[0] << "    \t"  << y[1] << "    \t"         
			<< y[2] << "    \t" << y[3]
			<< "    \t" << energy(y)     
			<< "    \t" << energy(y)-E0 
			<< "    \t" << h           
		        << "\n" ; 



Two ways you can do this... if you don't need to keep all the values returned by energy() as you vary y you can just test each new value to see if it's bigger or smaller than some max min and if it is then assign that new max min to max min. If you want to keep all the values assign each value to an element in an array and then loop through the array again with the same testing idea i laid out previously.


got it to work thanks

if (mn > energy(y)) {
     mn = energy(y);
    }
inserted in line 78(+ variable mn delaration at the start) finds the minimum

im sure i tried something like this, but i couldnt get it to work, my problem turned out to be pretty easy, thanks to everyone for the quick replies.

This post has been edited by araven7: 26 October 2010 - 06:12 PM

Was This Post Helpful? 0
  • +
  • -

#11 lfaivor21  Icon User is offline

  • New D.I.C Head

Reputation: 10
  • View blog
  • Posts: 29
  • Joined: 17-October 10

Re: Min/Max

Posted 26 October 2010 - 06:08 PM

Alright :)

This post has been edited by lfaivor21: 26 October 2010 - 06:13 PM

Was This Post Helpful? 0
  • +
  • -

#12 oscode  Icon User is offline

  • D.I.C Regular

Reputation: 109
  • View blog
  • Posts: 257
  • Joined: 24-October 10

Re: Min/Max

Posted 27 October 2010 - 03:19 AM

If you are using VS, Ctrl+A,Ctrl+K,Ctrl+F to add some better formatting to your code. Try to be consistent to make your code easier to read and maintain.
Was This Post Helpful? 0
  • +
  • -

Page 1 of 1