Particle Velocity in Fluid

Problem

Solution Steps and Example Problems

Analytical

Numerical

Figure 1: Roots of Equations f(v) = 0

Steps:

Substitute equation 1 into equation 4
The substitute result from step 1 into equation 3
Fill in all values provided in (b), such as gravity, diameter, etc.
Subtract" v" over to right side
Set "0" equal to f(v)

To calculate the root, observe where the function crosses the x-axis.

Sample Calculation

See below
Plug the "true" value into the combined equation to observe if the output is near 0.
Result: 0.0000185

Pseudocode Modified Secant Method

%% Set up and variable declaration%% Modified Secant Approach
%first estimatev0 = ((g)*(rho_s - rho_f)*(d^2))/((18)*(mu));v_sec(1) = v0;
DO while es < ea_sec(i) && (i<maxit) %Calculating function at xi %Calculating function at xi + delta*xi %Updates next value, or "xr". xr_old aka xi is effectively v_sec(i) here %True error %Approx error %Increment iend

Pseudocode Fixed-Point Method

%% Set up and variable declaration%% Modified Secant Approach
%first estimatev0 = ((g)*(rho_s - rho_f)*(d^2))/((18)*(mu));v_fp(1) = v0;
DOwhile (k<maxit) && (ea_fp(k) > es) v_fp(k+1) = sqrt((d_co)/(3*((24/... ((rho_f*v_fp(k)*d)/mu)) + (3/((sqrt((rho_f*v_fp(k)*d)/mu))))... + 0.34)*rho_f)); et_fp(k+1) = 100*(abs(true - v_fp(k+1))/true); ea_fp(k+1) = 100*(abs(v_fp(k+1)-v_fp(k))/v_fp(k+1)); k = k+1; end

Figure 2: Depicts Convergence and Velocity Equation Together

Other Important Results

Results

Looking at figure 2, we observe the root of the f(x) function coincides with the x-value of intersetion between the velocity equation and y=x. The original velocity equation (blue line) has a slope smaller than 1 near its intersection with y=x, which means the fixed-point method we used will converge towards the correct answer for any positive starting value.

The secant method required only four iterations to drop below the error requirement of 0.05%, whereas the fixed-point method took nine iterations.

We see that the Re value is above the 0.1 requirement to be considered laminar flow, meaning we must use the turbulant flow equations.

Error

Figure 3: True Percentage Error for Both Methods

Our error threshold "es" is below 0.05%. Looking at figure 3, we see that the secant method pulls ahead and completes the estimation of the root within four iterations. The fixed point method performs well in the initial iterations, but appears to approach 0 more slowly after four iterations. This error percentage is calculated by taking the difference between the true value observed in the analytical solution and the estimated value for the current iteration. It then divides this difference by the true value.

Figure 4: Approximate Percentage Error for Both Methods

Our error threshold "es" is below 0.05%. Looking at figure 4, we see that the secant method pulls ahead and completes the estimation of the root within four iterations. Interestingly, the fixed point method seems to outperform for the first iteration but subsequently takes more iterations to approach our error target.

Appendix I: MATLAB Code (expand to see code)

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%Author: Bryan Bennett%Date: 3/2/2020%Problem: Use the modified secant method and fixed-point iteration method to arrive%at a numerical understanding of the velocity of a particle in fluid. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% clear all; close all; %% Set up and variable declaration d = 0.0002 ; %Diameter (m)g = 9.80665; %Gravitational acceleration (m/s^2)rho_s = 7874; %Particle density (kg/m^3)rho_f = 1000; %Fluid density (kg/m^3)mu = 0.0014; %Dynamic viscosity of fluid (kg/m*s)delta = 10^-3; %Delta value for use in modified secant solutiones_error = 0.05; % (As a percentage) Error / stopping criteriondv = 0.001; %Step size in velocityfirst_v = 0; %Lower bound in velocity vectorlast_v = 0.5; %Upper bound in velocity vectorv = (first_v:dv:last_v); %Velocity vectormaxit = 1 + ((last_v - first_v)/dv); %Maximum number of iterations acrosss v vectori = 1; k = 1; %Initializing index parameterstrue = 59/847; %Represents "true value" of root %% Analytical Solution a_co = (72.*mu)./(v.*d);b_co = (9.*rho_f) ./ sqrt( (rho_f.*v.*d) ./ (mu) );c_co = 1.02*rho_f;d_co = 4.*g.*d.*(rho_s-rho_f); %This is f(v) = 0y_v = sqrt( (d_co) ./ ( a_co + b_co + c_co ) ) - v; Test_True = sqrt( (4*g*d*(rho_s-rho_f)) / ( ((72*mu)/(true*d)) + ((9*rho_f) / sqrt( (rho_f*true*d) / (mu) )) + (1.02*rho_f) ) ) - true;display(Test_True) %% Modified Secant Approach (answer to part b) %Initializing approx. and true error at 100%ea_sec(1) = 100; et_sec(1) = 100; %Equation 2 for estimating velocity under laminar conditions. Provides%first estimatev0 = ((g)*(rho_s - rho_f)*(d^2))/((18)*(mu)); % "xr = x0"v_sec(1) = v0; while es_error < ea_sec(i) && (i<maxit) %Calculating function at xi f_xi = sqrt( (d_co) / (3*(rho_f)*((24*mu)/(rho_f*v_sec(i)*d)+(3*sqrt(mu)*(sqrt(rho_f*v_sec(i)*d)))/... (rho_f*v_sec(i)*d)+0.34)))-v_sec(i); %Calculating function at xi + delta*xi f_xi_deltaxi = sqrt( (d_co) /(3*(rho_f)*((24*mu)/(rho_f*(v_sec(i)+(delta.*v_sec(i)))*d)+(3*sqrt(mu)*(sqrt(rho_f*(v_sec(i)+(delta.*v_sec(i)))*d)))/... (rho_f*(v_sec(i)+(delta.*v_sec(i)))*d)+0.34)))-(v_sec(i)+(delta.*v_sec(i))); %Updates next value, or "xr". xr_old aka xi is effectively v_sec(i) here v_sec(i+1) = v_sec(i) - ( (delta.*v_sec(i)*f_xi)/(f_xi_deltaxi-f_xi) ); %True error et_sec(i+1) = 100* (abs(true - v_sec(i+1))/true); %Approx error ea_sec(i+1) = 100* (abs(v_sec(i+1)-v_sec(i))/v_sec(i+1)); %Increment i i = i+1; end %Store final value / estimationsec_answer = v_sec(i); %% Reynolds Number (Re) and Drag Coefficient (Cd) (answer to part c) Re = (d*rho_f*v_sec(i))/(mu); %DimensionlessCd = (24/Re) + (3/(sqrt(Re))) + 0.34; %Dimensionless %% Fixed-Point Iteration Solution (answer to part d) ea_fp(1) = 100; et_fp(1) = 100; v_fp(1) = v0; while (k<maxit) && (ea_fp(k) > es_error) v_fp(k+1) = sqrt((d_co)/(3*((24/... ((rho_f*v_fp(k)*d)/mu)) + (3/((sqrt((rho_f*v_fp(k)*d)/mu))))... + 0.34)*rho_f)); et_fp(k+1) = 100*(abs(true - v_fp(k+1))/true); ea_fp(k+1) = 100*(abs(v_fp(k+1)-v_fp(k))/v_fp(k+1)); k = k+1; end fp_answer = v_fp(k); %% Plots (including answer to part e) figure(1)hold all;plot(v,y_v,'-.', 'LineWidth',3,'Markersize',6,'Color','g');plot(v,v,'-.', 'LineWidth',3,'Markersize',6,'Color','r');plot(v,y_v + v,'-.', 'LineWidth',3,'Markersize',6,'Color','b');xline(true,'r');xline(0); yline(0);xlabel('Function input (m/s)');ylabel('Function output (m/s)');title('Observing Convergence Alongside Root of Function');xlim([-.01 0.12]);ylim([-.01 0.1]);legend('f(v)','y = x','Original velocity equation', 'location', 'northwest')box on; grid on; i = (1:1:i);% figure (2)% hold all;% plot(i,ea_sec,'-.', 'LineWidth',3,'Markersize',6,'Color','g');% plot(i,et_sec,'-.', 'LineWidth',3,'Markersize',6,'Color','r');% xlabel('Number of Iterations');% ylabel('Percent Error {%}');% title('Percent Error for Modified Secant Method');% xlim([0 10]);% ylim([0 100]);% legend('Estimated Error (%)','True Error (%)','location', 'northeast')% box on; grid on; k = (1:1:k);% figure (3)% hold all;% plot(k,ea_fp,'-.', 'LineWidth',3,'Markersize',6,'Color','g');% plot(k,et_fp,'-.', 'LineWidth',3,'Markersize',6,'Color','r');% xlabel('Number of Iterations');% ylabel('Percent Error {%}');% title('Percent Error for Fixed Point Method');% xlim([0 10]);% ylim([0 100]);% legend('Estimated Error (%)','True Error (%)','location', 'northeast')% box on; grid on; figure (2)hold all;plot(i,ea_sec,'-.', 'LineWidth',3,'Markersize',6,'Color','g');plot(k,ea_fp,'-.', 'LineWidth',3,'Markersize',6,'Color','r');xlabel('Number of Iterations');ylabel('Percent Error {%}');title('Approximate Percent Error for both Methods');xlim([0 10]);ylim([0 100]);legend('Secant Method','Fixed Point','location', 'northeast')box on; grid on; figure (3)hold all;plot(i,et_sec,'-.', 'LineWidth',3,'Markersize',6,'Color','g');plot(k,et_fp,'-.', 'LineWidth',3,'Markersize',6,'Color','r');xlabel('Number of Iterations');ylabel('Percent Error {%}');title('True Percent Error for Both Methods');xlim([0 10]);ylim([0 100]);legend('Secant Method','Fixed Point','location', 'northeast')box on; grid on; figure(4)hold all;plot(v,y_v,'-o', 'LineWidth',3,'Markersize',2,'Color',[17 17 17]/255);xline(true,'r');xline(0); yline(0);xlabel('Function input (m/s)');ylabel('Function output (m/s)');title('Observing Root in Analytical Equation');xlim([-.01 0.12]);ylim([-.01 0.1]);legend('f(v)','location', 'northwest')box on; grid on;

Google Sites

Report abuse