Welcome to AE Resources
Converted document ROOT FINDING DEMONSTRATION

ROOT FINDING DEMONSTRATION

This program demonstrates the use of different root finding approaches to determine the local Mach number Ma when the pressure coefficient Cp, a and freestream Mach number M∞ are known. The nonlinear equation that described this relationship is:
(1) Cp, a = (2)/(γM2∞)[((1 + (γ − 1) ⁄ 2)M2∞)/(1 + (γ − 1) ⁄ 2)M2a))γ ⁄ (γ − 1) − 1] 
PROGRAM rootFinding
	IMPLICIT NONE
	REAL, PARAMETER :: pi=4.0*ATAN(1.0), gamma = 1.4
	INTEGER :: choice, nmax
	LOGICAL :: mystop = .TRUE.
	REAL :: M, Cp, Ma, ml, mr, tol, Cp_test, Cp_test2, slope
		
	!Input Mach number
	DO WHILE (mystop)
		WRITE(*,*) ’Enter the freestream Mach (M-Infinity) number’
		READ(*,*) M
		!Make sure Mach is reasonable
		IF ( M <= 0 .OR. M > 10) THEN
		  mystop = .TRUE.
		  WRITE(*,*) ’Too extreme of a Mach number.’
		  WRITE(*,*) ’Please choose a Mach number in  [0.1, 10].’
		ELSE
		  mystop = .FALSE.
		END IF
	END DO
	
	!Input Cp value
	DO
		WRITE(*,*) ’Enter coefficient of pressure (Cp) value’
		READ(*,*) Cp
		!Make sure this number is reasonable
		IF (Cp <= 1.0 .AND. Cp >= -10) THEN
			EXIT
		ELSE
			WRITE(*,*) ’Input Cp values between [-10,1]’
		END IF	
	END DO
	
	!Input numerical method choice
	DO
		WRITE(*,*) ’Enter <1> for Bisection Method or &
                    &<2> for Newton-Raphson Method’
		READ(*,*) choice
		!Makes sure user makes a legitimate choice
		IF (choice == 1 .OR. choice == 2) THEN
			EXIT
		ELSE
			WRITE(*,*) ’Please enter an appropriate choice (1 or 2)’
		END IF
	END DO
	
	nmax = 0  !This makes sure no infinite loops occur.
	SELECT CASE (choice)
		CASE (1)
			WRITE(*,*) ’’
			WRITE(*,*) ’Using the Bisection Method’
			!Set bounds
			ml = 0.0
			mr = 10.0
			DO
				!Set Ma at middle of bounds
				Ma = (ml + mr) / 2.0
				Cp_test = (2.0/(gamma*M**2)) * &
                     (( (1. + ((gamma - 1.)/2.)*(M**2.)) / &
             	 (1. + ((gamma - 1.)/2.)*(Ma**2.)) ) &
	          ** (gamma/(gamma - 1.)) - 1.) - Cp
				!Test tolerance
				tol = ABS(Cp_test)
				nmax = nmax + 1
				IF (tol <= .0001) THEN !Within tolerance of .1%
					EXIT
				ELSE IF (nmax >= 500) THEN !Stop infinite loop
					WRITE(*,*) ’Maximum counter limit reached (500).’
					WRITE(*,*) ’Error Analysis:’
					WRITE(*,*) ’     Either Mach or Cp may be too extreme, causing’
					WRITE(*,*) ’     Ma to never be found.  Please check input and retry.’
					STOP
				ELSE IF (Cp_test > 0) THEN
					!Resets left bound
					ml = Ma
				ELSE IF (Cp_test < 0) THEN
					!Resets right bound
					mr = Ma
				END IF
			END DO
Next Page →
Next Page →