Advertisement
Guest User

Untitled

a guest
Dec 8th, 2019
104
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
C++ 1.51 KB | None | 0 0
  1. #include <cmath>
  2. #include <iostream>
  3. #include <iomanip>
  4. using namespace std;
  5. double lambertW(double x, double prec=1e-8){
  6.     double w = 0;
  7.     for(int i=0;i<100;i++){
  8.         double wTimesExpW = w*exp(w);
  9.         double wPlusOneTimesExpW = (w+1)*exp(w);
  10.         w -= (wTimesExpW-x)/(wPlusOneTimesExpW-(w+2)*(wTimesExpW-x)/(2*w+2));
  11.         //cout <<w<<endl;
  12.         if (prec > abs((x-wTimesExpW)/wPlusOneTimesExpW)){
  13.             return w;
  14.  
  15.         }
  16.        
  17.            
  18.     /*  if (prec <= abs((x-wTimesExpW)/wPlusOneTimesExpW)){
  19.             cout <<"W(x) ne shoditsa x=" <<x<<endl;}*/
  20.     }
  21. return NAN;
  22. }
  23. double lambertW1(double x, double prec=1e-12){
  24.     double w=1.0;
  25.     for(int i=0;i<100;i++){
  26.         w=-log(-x/w);
  27.         //cout <<w<<endl;
  28.     }
  29.     return w;
  30. }
  31. double result(double x){
  32.         return  round(x*100000000)/100000000;
  33.    
  34. }
  35. int main()
  36. {
  37. double a,b;
  38. cin>>a>>b;
  39. if(a==0){
  40.     if(b<=0)
  41.     cout<<"no solution exists";
  42.     else
  43.     cout<<fixed<<setprecision(9)<<log(b);
  44.     return 0;
  45. }
  46. if(a==1&&b==1){
  47.     cout<<0;
  48.     return 0;
  49. }
  50. double x1=(-b/a)-lambertW(-1/(a*pow(M_E,b/a)));
  51. double x2=(-b/a)+lambertW1(-1/(a*pow(M_E,b/a)));
  52. x1=result(x1);
  53. x2=result(x2);
  54. if(x1==0)x1=abs(x1);
  55. if(x2==0)x2=abs(x2);
  56. if(isnan(x1)) {
  57.     if(isnan(x2)) {
  58.         cout<<"no solution exists";
  59.         return 0;
  60.     }else{
  61.         cout<<fixed<<setprecision(8)<<x2;
  62.         return 0;
  63.     }
  64. }else{
  65.     if(isnan(x2)) {
  66.         cout<<fixed<<setprecision(8)<<x1;
  67.         return 0;
  68.     }
  69. }
  70. if(x1<x2)
  71.     cout<<fixed<<setprecision(8)<<x1<<" "<<fixed<<setprecision(8)<<x2;
  72. else
  73.     cout<<fixed<<setprecision(8)<<x2<<" "<<fixed<<setprecision(8)<<x1;
  74. return 0;
  75. }
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement