# Python Program to Implement Steady State Shock Equations Assignment Solution

July 03, 2024
Martin Jonas
🇦🇺 Australia
Python
Dr. Martin Jonas, PhD in Computer Science from Southern Cross University, Australia. With 4 years of experience in Python assignments, I offer expert guidance and support to help you excel in your programming projects.
Key Topics
• Instructions
• Objective
• Requirements and Specifications
Tip of the day
News

## Instructions

### Objective

Write a python homework to implement steady state shock equations.

## Requirements and Specifications

Source Code

import numpy as np import matplotlib.pyplot as plt M = 6.0 gam = 5./3. # ========================================================================== # ========================================================================== def ddx(x, F): Ns = F.size - 1 dFdx = np.zeros(Ns + 1, dtype='double') dFdx[0] = (F[1] - F[0]) / (x[1] - x[0]) for i in range(1, Ns): dFdx[i] = (F[i + 1] - F[i - 1]) / (x[i + 1] - x[i - 1]) dFdx[Ns] = (F[Ns] - F[Ns - 1]) / (x[Ns] - x[Ns - 1]) return dFdx # ========================================================================== # ========================================================================== if __name__ == '__main__': electron = np.loadtxt('./electron_quantites.txt') ion = np.loadtxt('./ion_quantites.txt') x = electron[:,0] # positions rho = ion[:,1] # this is rho/rho0. ke = electron[:,2] Te = electron[:,1] Ti = ion[:,2] ki = ion[:,3] eta = ion[:,4] u = ion[:,5] p = ion[:,6] E = ion[:,7] # First, calculate du/dx dudx = ddx(x, u) # Calculate dTi/dx dTidx = ddx(x, Ti) # Calculate dTe/dx dTedx = ddx(x, Te) # From equation (1b), calculate p0u0^2 + p0 eq1b = p + gam*M**2.*rho*u**2. - eta*dudx eq1b_right = gam*M**2 + p[0] # Equation (1c) eq1c = E*u + p*u - ki*dTidx - ke*dTedx - eta*u*dudx eq1c_right = E[0]/rho[0] + p[0] # Eq. 1a eq1a = rho*u #========================================================================== #========================================================================== # create plot fig, axes = plt.subplots(nrows = 2, ncols = 3, figsize=(7,7)) # Plot dTidx axes[0,0].plot(x, rho) axes[0,0].set_title(r'$\rho$') axes[0,0].set_xlabel('Position') axes[0,0].set_ylabel('Plasma mass density') axes[0,0].grid(True) axes[0,1].plot(x, ke, label = 'ke') axes[0,1].set_title(r'$k_{e}$') axes[0,1].set_xlabel('Position') axes[0,1].set_ylabel('Electron Conductivity') axes[0, 1].grid(True) axes[0,2].plot(x, Te, label = 'Te') axes[0,2].set_title(r'$T_{e}$') axes[0,2].set_xlabel('Position') axes[0,2].set_ylabel('Electron Temperature') axes[0, 2].grid(True) axes[1,0].plot(x, Ti, label = 'Ti') axes[1,0].set_title(r'$T_{i}$') axes[1,0].set_xlabel('Position') axes[1,0].set_ylabel('Ion Conductivity') axes[1, 0].grid(True) axes[1,1].plot(x, ki, label = 'ki') axes[1,1].set_title(r'$k_{i}$') axes[1, 1].set_xlabel('Position') axes[1, 1].set_ylabel('Ion conductivity') axes[1, 1].grid(True) axes[1,2].plot(x, eta, label = 'eta') axes[1,2].set_title(r'$\eta$') axes[1, 2].set_xlabel('Position') axes[1, 2].set_ylabel('Ion viscosity') axes[1, 2].grid(True) plt.show() fig, axes = plt.subplots(nrows = 2, ncols = 3, figsize=(7,7)) axes[0,0].plot(x, u, label = 'u') axes[0,0].set_title(r'$u$') axes[0, 0].set_xlabel('Position') axes[0, 0].set_ylabel('Bulk fluid velocity') axes[0, 0].grid(True) axes[0,1].plot(x, p, label = 'p') axes[0,1].set_title(r'$p$') axes[0, 1].set_xlabel('Position') axes[0, 1].set_ylabel('Plasma pressure') axes[0, 1].grid(True) axes[0,2].plot(x, E, label = 'E') axes[0,2].set_title(r'$E$') axes[0, 2].set_xlabel('Position') axes[0, 2].set_ylabel('Total Plasma Energy') axes[0, 2].grid(True) axes[1,0].plot(x, dudx, label = 'dudx') axes[1,0].set_title(r'$\frac{\partial u}{\partial x}$') axes[1, 0].set_xlabel('Position') axes[1, 0].set_ylabel('Rate of change of bulk velocity') axes[1, 0].grid(True) axes[1,1].plot(x, dTidx, label = 'dTidx') axes[1,1].set_title(r'$\frac{\partial T_{i}}{\partial x}$') axes[1, 1].set_xlabel('Position') axes[1, 1].set_ylabel('Rate of change of Ion Temp.') axes[1, 1].grid(True) axes[1,2].plot(x, dTedx, label = 'dTedx') axes[1,2].set_title(r'$\frac{\partial T_{e}}{\partial x}$') axes[1, 2].set_xlabel('Position') axes[1, 2].set_ylabel('Rate of change of Electron Temp.') axes[1, 2].grid(True) plt.show() fig = plt.figure(figsize=(7, 5), dpi=120) plot = plt.plot(x, (eq1b-eq1b_right)/eq1b_right*100., 'royalblue', linestyle='--', label='Eq. (1b) percent error') plot = plt.plot(x, (eq1c-eq1c_right)/eq1c_right*100., 'tomato', linestyle='-', label='Eq. (1c) percent error') plot = plt.plot(x, (eq1a - 1.0)/1.0 *100, 'purple', linestyle='-', label = 'Eq. (1a) percent error') plt.ylabel(r'Percent error [%]') plt.xlabel(r'Position $\hat{x}$') #plt.xlim([25.0, 32.5]) plt.xlim([min(x), max(x)]) legend = plt.legend(loc='best', shadow=False, fontsize='small') #plt.savefig('equation_error.eps', format='eps', dpi=1000) plt.grid(True) plt.show() plt.close() #========================================================================== #==========================================================================

## Similar Samples

Explore our diverse range of programming homework samples across Java, Python, C++, and more. Each sample exemplifies our proficiency in solving complex coding challenges with precision and clarity. Whether you're seeking help with algorithms, data structures, or software development projects, our samples showcase our commitment to delivering high-quality solutions tailored to your academic needs.