Advertisement
Guest User

Ternary Search

a guest
Apr 25th, 2019
115
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. #include <bits/stdc++.h>
  2. using namespace std ;
  3.  
  4. double M_a , M_b ,  V_a , V_eq , K_a , K_w , V_m , V_p ;
  5.  
  6. double pH(double V_b)
  7. {
  8.     double beta , b , c , x , E ;
  9.     x = (V_eq+(K_a/M_b)*(V_a+V_eq)) ;
  10.  
  11.     beta = (K_w*K_w*(V_a+V_eq)*(V_a+V_eq))/(2*M_b*M_b*x*x) ;
  12.     beta *= (1 + sqrt(1 + (4*K_a*M_b*x)/(K_w*(V_a+V_eq)))) ;
  13.     beta += (K_a*K_w*(V_a+V_eq))/(M_b*x) ;
  14.     beta = 1/sqrt(beta) ;
  15.  
  16.     x = V_b + (K_a*(V_a+V_b))/M_b ;
  17.     c = -(K_a*K_w*(V_a+V_b))/(M_b*x) ;
  18.     b = (K_a*(V_b - (M_a/M_b)*V_a)) - (K_w/M_b)*(V_a+V_b) ;
  19.     b /= x ;
  20.  
  21.     x = -(b/2) + sqrt(((b*b)/4)-c) ;
  22.     E = ((beta*beta*x*x)-1)/(2*beta*x) ;
  23.  
  24.     double ret = log10(beta) - (asinh(E)*log10(exp(1))) ;
  25.    
  26.     return ret ;
  27. }
  28.  
  29. double err(double u){
  30.     M_a = u ;
  31.     V_eq = (M_a*V_a)/M_b ;
  32.     double mm = (3.1+4.4)/2 ;
  33.     double mp = (8.2+10)/2 ;
  34.     double ret = abs(mm - pH(V_m))+abs(mp-pH(V_p)) ;
  35.     return ret ;
  36. }
  37.  
  38. int main(int argc, char const *argv[])
  39. {
  40.     K_w = 1e-14 ;
  41.  
  42.     printf("Enter the volume of the weak acid: ") ;
  43.     scanf("%lf", &V_a) ;
  44.     printf("Enter the equilibrium constant of the weak acid: ") ;
  45.     scanf("%lf", &K_a) ;
  46.     printf("Enter the concentration of the strong base: ") ;
  47.     scanf("%lf", &M_b) ;
  48.     printf("Enter the volume of the strong base in case of Methyl Orange: ") ;
  49.     scanf("%lf", &V_m) ;
  50.     printf("Enter the volume of the strong base in case of Phenolphthalein: ") ;
  51.     scanf("%lf", &V_p) ;   
  52.  
  53.     //ternary search
  54.     double l = 0 , r = 1 ;
  55.    
  56.     while(r-l > 1e-6){
  57.         double l1 = (l*2+r)/3 ;
  58.         double l2 = (l+2*r)/3 ;
  59.        
  60.        if(err(l1) < err(l2))r = l2 ;
  61.        else l = l1 ;
  62.        
  63.     }
  64.  
  65.     printf("%lf\n" , l) ;
  66.  
  67.     return 0;
  68. }
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement