|
[Sponsors] |
Falkner skan solution with different beta value |
|
LinkBack | Thread Tools | Search this Thread | Display Modes |
October 10, 2015, 17:16 |
Falkner skan solution with different beta value
|
#1 |
New Member
Ankur Bhatnagar
Join Date: Nov 2013
Posts: 19
Rep Power: 12 |
I am trying to write falkner skan code with varying value of beta(2*m/m+1), where m tell us about velocity profile. My code looks like...
#include<iostream> #include<cstdlib> #include<cmath> #include<iomanip> #include<fstream> using namespace std; double function_f1(double z,double z1,double z2,double eta); double function_f2(double z,double z1,double z2 ,double eta); double function_f3( double z,double z1,double z2,double eta,double beta); int main(){ double beta,deta,res=10.0; double *z,*z1,*z2,*z2new; double *eta; int e,max_e=7001,iter; double k[4],m[4],n[4]; double A,A_improved; eta=new double [max_e+1]; z=new double [max_e+1]; z1=new double [max_e+1]; z2=new double [max_e+1]; z2new=new double [max_e+1]; //initial values eta[0]=0.0L; deta= 0.001L; z[0]=0.0L; z1[0]=0.0L; z2[0]=0.322L; A=z2[0]; double *diff,*err,*q; diff=new double [50]; err =new double [50]; q =new double [50]; beta = -0.0 ; iter=0; for (int iter=0;iter<12;iter++) { for (e=0;e<max_e;e++) { k[0]= function_f1(z[e],z1[e],z2[e],eta[e]); m[0] = function_f2(z[e],z1[e], z2[e],eta[e]); n[0] = function_f3(z[e],z1[e],z2[e],eta[e], beta); k[1]= function_f1(z[e]+0.5*k[0]*deta,z1[e]+0.5*m[0]*deta,z2[e]+0.5*n[0]*deta,eta[e]+0.50*deta); m[1] = function_f2(z[e]+0.5*k[0]*deta,z1[e]+0.5*m[0]*deta,z2[e]+0.5*n[0]*deta,eta[e]+0.50*deta); n[1] = function_f3(z[e]+0.5*k[0]*deta,z1[e]+0.5*m[0]*deta,z2[e]+0.5*n[0]*deta,eta[e]+0.50*deta,beta); k[2]= function_f1(z[e]+0.5*k[1]*deta,z1[e]+0.5*m[1]*deta,z2[e]+0.5*n[1]*deta,eta[e]+0.50*deta); m[2] = function_f2(z[e]+0.5*k[1]*deta,z1[e]+0.5*m[1]*deta,z2[e]+0.5*n[1]*deta,eta[e]+0.50*deta); n[2] = function_f3(z[e]+0.5*k[1]*deta,z1[e]+0.5*m[1]*deta,z2[e]+0.5*n[1]*deta,eta[e]+0.50*deta,beta); k[3]= function_f1(z[e]+k[2]*deta,z1[e]+m[2]*deta,z2[e]+n[2]*deta,eta[e]+0.50*deta); m[3] = function_f2(z[e]+k[2]*deta,z1[e]+m[2]*deta,z2[e]+n[2]*deta,eta[e]+0.50*deta); n[3] = function_f3(z[e]+k[2]*deta,z1[e]+m[2]*deta,z2[e]+n[2]*deta,eta[e]+0.50*deta,beta); z[e+1] = z[e]+ (deta/6.0)* ( k[0] + 2.0*k[1] +2.0*k[2] + k[3] ) ; z1[e+1] = z1[e]+ (deta/6.0)* ( m[0] + 2.0*m[1] +2.0*m[2] + m[3] ); z2[e+1] = z2[e]+ (deta/6.0)* ( n[0] + 2.0*n[1] +2.0*n[2] + n[3] ); eta[e+1] =eta[e]+ deta; //diff[iter]= z1[1088]; //random point near 1.0 choosen from the curve } ofstream outfile; outfile.open("falknerbeta0.0.dat"); for (e=0;e<max_e;e++) { cout<< eta[e] <<" "<<z[e]<<" "<<z1[e]<<" "<<z2[e]<<endl; outfile<< eta[e] <<" "<<z[e]<<" "<<z1[e]<<" "<<z2[e]<<endl; } outfile.close(); /*//newton raphson q[iter]= A; err[iter] = (1.0L - z1[1088]); if (iter == 1) {A=A+ 0.05;} else if (iter == 3) if (err[2] <err[1]) A_improved=q[2]-((q[3]-q[2])*err[2]/(err[3]-err[2])); else {A_improved=q[1]-((q[3]-q[1])*err[1]/(err[3]-err[1]));} else { A_improved=q[iter-1]-((q[iter]-q[iter-1])*err[iter-1]/(err[iter]-err[iter-1]));} {q[iter]=q[iter-1]-((q[iter]-q[iter-1])*err[iter-1]/(err[iter]-err[iter-1]));} if (fabs(err[iter]> 0.0001)){ z2[0]= q[iter];} else {exit(0);}*/ z2[0]=z2[0]+0.001; } delete z2new; delete z2; delete z1; delete z; return 0; } double function_f1(double z,double z1,double z2,double eta){ double f1; f1= z1; return f1; } double function_f2(double z,double z1,double z2,double eta){ double f2; f2= z2; return f2; } double function_f3( double z,double z1,double z2,double eta,double beta){ double f3; f3= -0.5*z*z2 - beta * (1.0L-z1*z1); return f3; } i am using RK4 with values try to refine using newton raphson... My question is 1.if i am doing it correctly. Is there any better way to write this code in c++. 2. for blasius profile it is giving correct result but become unstable at higher beta value(beta =1.0, 0.5 ..) .why?(i dont know) 3. the commented part is newton raphson ...i tried to refine my rk4 solution with NR but could not implement it correctly...any one please try thanks ! |
|
|
|
Similar Threads | ||||
Thread | Thread Starter | Forum | Replies | Last Post |
grid dependancy | gueynard a. | Main CFD Forum | 19 | June 27, 2014 21:22 |
Why 3D solid-pore geometry showing diverged solution? | Sargam05 | OpenFOAM | 0 | December 3, 2012 15:45 |
Analytic solution for 2D steady Euler equations | jojo81 | Main CFD Forum | 0 | October 15, 2012 12:05 |
Determining alpha and beta for porous baffle | Liaqat Khan | Siemens | 1 | October 27, 2000 04:44 |
Wall functions | Abhijit Tilak | Main CFD Forum | 6 | February 5, 1999 01:16 |