hello admin, please help check my C++ code. it complied but didnt return 0, so n output. The codes are below; //HOWARTH TRANSFORMATION OF BLASIUS BOUNDARY LAYER PRROBLEM // USING 100STEPS i.e t=0.05, //ASSUMING INFINITY=5
#include<iostream> #include "stdafx.h" #include <stdlib.h> #include <stdio.h> #include <math.h> #include "solvers.h" #include<cmath> #include<iomanip>
using namespace std; using namespace alglib;
int main(int argc, char **argv) { int m,n,i,j,k,M;// sum={0}; double A[300][1], B[300][300], C[300][1], D[300][1], E[300][1], F[300][1], s[300]={0}, f[300], h[300], T[300] ,t,a,b,w,z,y,x, p;// l[300][300]={0},u[300][300]={0},z[300]={0}; cout<<"Aim: To solve Howarth Transformation Boundary Layer problem using nonlinear solver (Newton Raphson method) for: \n"<<" f'''+ f/(h^(-1/3)) f''= 0 \n"<<" h''+ (0.7fh'/h^(-1/3))- 0.28M^(2)*f''^(2)=0 \n"<<endl; cout<<"Enter number of unknown variables "<<endl; cin>>m; cout<<"Enter number of non-linear equations "<<endl; cin>>n; cout<<"Enter initial guess value for f "<<endl; cin>>a; cout<<"Enter initial guess value of for h "<<endl; cin>>b; for(i=1;i<=99;i++) { f[i]= a; h[i]= b; }
/* 198 non-linear functions with 198 unknowns Fn(f1,f2,f3,f4,h1,h2,h3,h4) */ A[0][0]= -2*f[1] - 2*f[2] + f[3] + (0.1*f[1]*f[2])/pow(h[1],-1/3) - ((0.2*pow(f[1],2))/pow(h[1],-1/3)); A[1][0]= 2*f[1] - 2*f[3] + f[4] + (0.1*f[1]*f[2])/pow(h[2],-1/3) + (0.1*f[2]*f[3])/pow(h[2],-1/3) - (0.2*pow(f[2],2))/pow(h[2],-1/3); for(i=3;i<=97;i++) { A[i-1][0]= f[i+2]-2*f[i+1]+2*f[i-1]-2*f[i-2] + 0.1*f[i]*f[i+1]/pow(h[i],-1/3) + 0.1*f[i]*f[i-1]/pow(h[i],-1/3) - (0.2*pow(f[i],2))/pow(h[i],-1/3); } A[97][0]= 0.1-f[99]+ 2*f[97]-2*f[96] + (0.1*f[98]*f[99])/pow(h[98],-1/3) +(0.1*f[98]*f[97])/pow(h[98],-1/3) - ((0.2*pow(f[98],2))/pow(h[98],-1/3)); A[98][0]= -f[99]+ 2*f[98]-2*f[97] + (0.005*f[99])/pow(h[99],-1/3) + (0.1*f[99]*f[98])/pow(h[99],-1/3) - (0.1*pow(f[99],2))/pow(h[99],-1/3) + 0.05; A[99][0]= pow(M,2)*(-2.24*pow(f[1],2) - 0.56*pow(f[2],2) + 2.24*f[1]*f[2] )- 0.01*h[1] + 0.005*h[2] + (0.0000875*f[1]*h[2])/pow(h[1],-1/3) - (0.000175*f[1])/pow(h[1],-1/3) + 0.01; for(i=2; i<=98; i++) { A[98+i][0]= pow(M,2)*(-2.24*pow(f[i],2) - 0.56*pow(f[i-1],2) - 0.56*pow(f[i+1],2) + 2.24*f[i]*f[i-1] - 1.12*f[i-1]*f[i+1] + 2.24*f[i]*f[i+1] )- 0.01*h[i] + 0.005*h[i-1] + 0.005*h[i+1] + (0.0000875*f[i]*h[i+1])/pow(h[i],-1/3) - (0.0000875*f[i]*h[i-1])/pow(h[i],-1/3); } A[197][0]= pow(M,2)*(-0.0014 -0.0056*f[98] + 0.0056*f[99] - 0.56*pow(f[98],2) - 0.56*pow(f[99],2) + 1.12*f[98]*f[99]) + 0.005*h[98] - 0.01*h[99] + (0.0000875*f[99])/pow(h[99],-1/3) - (0.0000875*f[99]*h[98])/pow(h[99],-1/3) + 0.005 ;
/* Jacobian Matrix elements of Matrix B */ B[0][0]= -2+(0.1*f[2])/pow(h[1],-1/3)-(0.4*f[1])/pow(h[1],-1/3); B[0][1]= -2+(0.1*f[1])/pow(h[1],-1/3); B[0][2]= 1; B[0][3]= 0; B[0][99]= (0.1*f[1]*f[2]*pow(h[1],-2/3)/3)-(0.2*pow(f[1],2)*pow(h[1],-2/3)/3); B[0][100]= 0; B[0][101]= 0; B[0][102]= 0; B[1][0]= 2+(0.1*f[2])/pow(h[2],-1/3); B[1][1]= (0.1*f[1])/pow(h[2],-1/3)+ (0.1*f[3])/pow(h[2],-1/3)-(0.4*f[2])/pow(h[2],-1/3); B[1][2]= -2+ (0.1*f[2])/pow(h[2],-1/3) ; B[1][3]= 1; B[1][99]= 0; B[1][100]= (0.1*f[1]*f[2]*pow(h[2],-2/3)/3) + (0.1*f[2]*f[3]*pow(h[2],-2/3)/3) - (0.2*pow(f[2],2)*pow(h[2],-2/3)/3); B[1][101]= 0; B[1][102]= 0; for(i=3; i<=97; i++) { B[i-1][i-1]= 0.1*f[i+1]/pow(h[i],-1/3) - 0.4*f[i]/pow(h[i],-1/3) + 0.1*f[i-1]/pow(h[i],-1/3); B[i-1][i-2]= 2 + 0.1*f[i]/pow(h[i],-1/3); B[i-1][i+1] = 1 ; B[i-1][i]= -2 + 0.1*f[i]/pow(h[i],-1/3); B[i-1][i-3] = -2 ; B[i-1][98+i]= (0.1*f[i]*f[i+1]*(h[i],-2/3)/3) - (0.2*(f[i],2)*(h[i],-2/3)/3) + (0.1*f[i]*f[i-1]*(h[i],-2/3)/3); }
B[97][95]= -2; B[97][96]= 2+(0.1*f[98])/pow(h[98],-1/3); B[97][97]= (0.1*f[99])/pow(h[98],-1/3) + (0.1*f[97])/pow(h[98],-1/3)- (0.4*f[98])/pow(h[98],-1/3); B[97][98]= -1 + (0.1*f[98])/pow(h[98],-1/3); B[97][196]= (0.1*f[98]*f[99]*pow(h[98],-2/3)/3) + (0.1*f[98]*f[97]*pow(h[98],-2/3)/3)-(0.2*pow(f[98],2)*pow(h[98],-2/3)/3); B[98][96]= -2; B[98][97]= 2+(0.1*f[99])/pow(h[99],-1/3); B[98][98]= -1+ (0.005)/pow(h[99],-1/3)+ (0.1*f[98])/pow(h[99],-1/3)-(0.2*f[99])/pow(h[99],-1/3); B[98][197] = (0.1*f[99]*f[98]*pow(h[99],-2/3)/3) + (0.01*f[99]*pow(h[99],-2/3)/3)-(0.1*pow(f[99],2)*pow(h[99],-2/3)/3); B[99][0]= pow(M,2)*(-4.48*f[1] + 2.24*f[2]) + (0.0000875*h[2])/pow(h[1],-1/3) - (0.000175)/pow(h[1],-1/3) ; B[99][1]= pow(M,2)*(2.24*f[1] - 1.12*f[2]); B[99][2]= 0; B[99][3]= 0; B[99][99]= -0.01 + (0.0000875*f[1]*h[2]*pow(h[1],-2/3)/3) - ( 0.000175*f[1]*pow(h[1],-2/3)/3); B[99][100]= 0.005 + (0.0000875*f[1])/pow(h[1],-1/3); B[99][101]= 0; B[99][102]= 0; for(i=2; i<=98; i++) { B[98+i][i-1]= pow(M,2)*(2.24*f[i+1] - 4.48*f[i] + 2.24*f[i-1] ) + (0.0000875*h[i+1])/pow(h[i],-1/3) - (0.0000875*h[i-1])/pow(h[i],-1/3) ; B[98+i][i-2]= pow(M,2)*(-1.12*f[i+1] + 2.24*f[i] - 1.12*f[i-1]); B[98+i][i]= pow(M,2)*(-1.12*f[i+1] + 2.24*f[i] - 1.12*f[i-1] ); B[98+i][98+i]= -0.01 + (0.0000875*f[i]*h[i+1]*pow(h[i],-2/3)/3 ) - (0.0000875*f[i]*h[i-1]*pow(h[i],-2/3)/3 ); B[98+i][97+i]= 0.005 - (0.0000875*f[i])/pow(h[i],-1/3) ; B[98+i][99+i]= 0.005 + (0.0000875*f[i])/pow(h[i],-1/3) ; } B[197][97]= pow(M,2)*(-0.0056 - 1.12*f[98] + 1.12*f[99] ); B[197][98]= pow(M,2)*(0.0056 - 1.12*f[99] + 1.12*f[98] ) + (0.0000875/pow(h[99],-1/3) - (0.0000875*h[98])/pow(h[99],-1/3)); B[197][196]= 0.005 - (0.0000875*f[99])/pow(h[99],-1/3) ; B[197][197]= -0.01 + (0.0000875*f[99]*pow(h[99],-2/3)/3) - ( 0.0000875*f[99]*h[98]*pow(h[99],-2/3)/3);
//for(i=0;i<198;i++) // for(j=0;j<198;j++) // { // ae_int_t w = 198; // sparsematrix B[i][j]; // real_1d_array z= "A[i][j]"; // real_1d_array s[i]; // sparsesolverreport rep; // bool isuppertriangle = false; // sparsesolvesks(B[i][j], w, isuppertriangle, z, rep, s[i]); //} // printf("%s\n", s[i]); //
for(i=0;i<198;i++) for(j=0;j<198;j++) { real_2d_array y= "B[i][j]"; ae_int_t n = 198; bool isupper = false; real_2d_array z = "A[i][j]"; ae_int_t m = 198; ae_int_t info; densesolverreport rep; real_2d_array s[i]; smp_spdmatrixsolvem(y, n, isupper, z, m, info, rep, s[i]); } printf("%s\n", s[i]); return 0; }
|