More sample programs are coming up: Delaunay triangulation code, a panel code for an airfoil, 2D unstructured Navier-Stokes code, etc.
[Note: As the name, Katate Masatsuka, implies, I write only when I find time.]

"I do like CFD, VOL.1, Second Edition" is now available in both printed and PDF versions. Download PDF (FREE) at cfd-boook page.

CFD seminar at National Institute of Aerospace is now broadcast online!
Check out the website for seminar videos and files - NIA CFD Seminar

Favorite Codes and Subroutines      

---

Educationally-Designed Unstructured 2D (EDU2D) Download all files: edu2d_twod2threed.zip This program reads a 2D tria/quqad/mixed grid, and generates a 3D grid by extending/rotating the 2D grid to the third dimension. A hexahedral grid will be generated for a pure quadrilateral grid, and a prismatic or tetrahedral grid will be generated for a pure triangular grid as below. Below, a flat-plate grid is used as an example. It generates a .ugrtid file (for a solver such as FUN3D), .su2 file for SU2, and tecplot and vtk files for viewing the grid.

3D conventional and hyperbolic Poisson solver:   edu3d_hyperbolic_poisson_v1.tar.gz" This is a 3D Poisson solver code having two options: (1) conventional solver and (2) hyperbolic solver. The hyperbolic solver is faster for finer grids than the conventional solver beause the numerical stiffness caused by second-derivatives in the Poisson equation is completely elimianted in the hyperbolic solver; and it is also more accurate on unstructured grids, especially in the solution derivatives (one order higher accurate than conventional methods). See JCP2020 for the details of the conventional and hyperbolic solvers. Example cases are included, where linear and nonlinear Poisson equations are solved in a cube domain and also in a domain bounded by two co-centroid spheres, using tetrahedral and hexahedral grids as shown in figures below.

3D Grid Generation for Sphere:   edu3d_sphere.tar.gz A 3D grid is generated for a sphere and written in the UGRID format. A tetrahedarl, prismatic, or mixed grid can be generated with a specified y-plus value. Learn how a 3D unstructured grid is generated over a sphere.

Educationally-Designed Unstructured 2D (EDU2D) Exact vortex solution to the 2D unsteady Euler equations:     edu2d_inviscid_vortex_solution_verification.f90 This program computes the inviscid vortex solution and substitute it into the 2D unsteady Euler equations and show all equations are zero (satisfied).

Educationally-Designed Unstructured 2D (EDU2D) EDU2D Cylinder-Wake II: Grid generation for cylinder:   Download all files: 'edu2d_cylinder_wake_version02.tar.gz This is another grid generation code for a domain around a circular cylinder. It generates multiblock-like grid with better controls than the other code. This package contains the twod2threed code and generates a 3D grid by extending the 2D (x,z) grid in the y-direction. For the grid format, see   edu2d_grid_format.pdf

Educationally-Designed Unstructured 3D (EDU3D) EDU3D_MMS_Euler   edu3d_mms_euler_v1.f90 This code shows how to compute the source terms in the method of manufactured solutions (MMS) for the 3D Euler equations. Any function can be made an exact solution to the 3D Euler equations with suitable source terms. That is, any function v(x,y,z) is an exact solution to the following equation: MMS has been widely used for verifying the order of accuracy of a CFD code. Grab this code, and learn how it can be implemented. I always use MMS to verify the implementation. It is a very useful tool for debugging.

Educationally-Designed Unstructured 3D (EDU3D) EDU3D Euler:   Download all files: edu3d_euler.tar.gz This is an implicit unstructured-grid Euler solver (f90, serial, tetra only), written for an educational purpose. Read the source code and learn how a node-centered edge-based unstructured-grid solver is written. Example cases are included: MMS test, cube, hemisphere cylinder, and ONERA M6 wing. For the grid format, see   edu2d_grid_format.pdf

Heapsort program   heapsort_v1.f90 I have used the heapsort algorithm in a hierarchical agglomeration scheme developed and implemented in NASA's FUN3D code. This is just an example program of the heapsort. In FUN3D, I used the heapsort to order the list of control-volumes to be agglomerated based on the number of non-agglomerated neighbor volumes, so that agglomeration can be performed from the volumes with the least # of neighbors, which helps keep the front as convex as possible and avoid creating holes. The front list changes every time a volume is agglomerated with neighbors and it needs to be re-ordered with new front members (the neighbors of the newly agglomerated volumes). It was so fast that the agglomeration can be peformed quickly inside the solver every time it reads a grid.

EDU2D-CCFV-EULER-EXPLCT: An explicit Euler code: Download all files: edu2d-ccfv-euler-explct_2019.zip This is an explicit unstructured-grid Euler code, written for an educaitonal purpose. Read the source code and learn how a cell-centered FV unstructured-grid solver is written. Examples cases are included: transonic airfoil, subsonic cylinder, unsteady vortex, numerical truncation error computation. For the grid format, see   edu2d_grid_format.pdf

EDU2D-Rectangular-GRID: Grid generation for a rectangular domain:   edu2d_rectangular_grid_v1.f90 This code generates a triangular or quadrilateral grid for a 2D rectangular domain with 4 or 5 boundary parts as shown below (Top: 4 boundary parts, Bottom: 5 boundary parts). The 5-boundary-part grid can be used for a shock diffraction problem or a forward-facing step problem, for example. It generates a custom 2D grid file (used in EDU2D codes), tecplot, SU2, and vtk files.

EDU2D-VISCOUS-BOX-GRID: Grid generation for a viscous-wall domain:   edu2d_viscous_box_grid_v0.f90 This code generates a triangular or quadrilateral grid for a 2D rectangular domain with 4 viscous boundaries as shown below, suitable for natural convection (buoyancy flow) problems. It generates a custom 2D grid file (used in EDU2D codes), tecplot, SU2, and vtk files.

EDU2D-FP-GRID: Grid generation for flat plate:   edu2d_fp_grid_v1.f90   input.nml This code generates a grid (tria,quad,mixed) for a 2D flat plate case, where a flat plate is placed on the bottom right half, (x,y)=(0,0)-(xmax,0). It generates a custom 2D grid file (used in EDU2D codes), tecplot, SU2, and vtk files.

EDU2D-HALF-CYLINDER-GRID: Grid generation for half cylinder:   edu2d_half_cylinder_grid_v0.f90

EDU2D-BUMP-GRID: Grid generation for bump:   edu2d_bump_irregular_grid_v1.f90

Flow around a Karman-Trefftz airfoil (generate grids and exact solutions):   edu2d_vkt_airfoil_v4.f90 ,   vkt_airfoil_v1_display.m Karman-Trefftz airfoil is an intriguing airfoil for which a complete set of exact solutions can be computed. This solution has been used by many people to verify the accuracy of their inviscid code (incompressible limit). If you have not, use this code now to generate a (quadrilateral or triangular) grid, run your code on it, compute the error, and verify the accuracy of your code. Me? I have used it for my third-order multigrid Cauchy-Riemann solver (IJNMF 2004). Updated (10-13-10). Formula corrected. It works now for cambered airfoils.

Ringleb's Flow (generate grids and exact solutions):   ringleb_v1.f90 Ringleb's flow is a famous exact solution of the compressible Euler equations with a smooth transition from subsonic to supersonic without any shock waves. This solution has been used actually by many people to verify the accuracy of their Euler code. If you have not, use this code now to generate a (quadrilateral or triangular) grid, run your Euler code on it, compute the error, and verify the accuracy of your Euler code. Minor bugs fixed (12-29-10).

Flow over Rankine's half body (generate grids and exact solutions):   edu2d_rankine_half_body_grid_v1.f90 A flow over Rankine's half body is an interesting exact solution to the potential flow equations: Cauchy-Riemann systems, Lapace equations, the potential equation, Euler equations. You can download, compile, and run this code to generate a (quadrilateral or triangular) grid, with the exact solution computed at each grid point (stored in Tecplot files). At some point, I'd like to use this solution for verification. Shown below is the exact x-velocity contours on an irregular triangular grid .

EDU2D-Cylinder-Wake: Grid generation for cylinder:   edu2d_cylinder_wake_v4.f90   input.nml [To compile: gfortran edu2d_cylinder_wake_v4.f90] This code generates a grid (tria,quad,mixed) over a cylinder with a higher-resolution in the wake region. It is suitable for unsteady cylinder simulations. It generates a custom 2D grid file (used in EDU2D codes), tecplot, SU2, and vtk files.

Educationally-Designed Unstructured 2D (EDU2D) Special triangular grid generator: EDU2D-50yen-Tria-Grids   edu2d_50yen_tria_grids_v2.f90   edu2d_50yen_tria_grids_tecplot_v0.lay This program can be useful for generating triangular grids for use in CFD or just learning how to manipulate unstructured-grid (invert, stretch, desired grid-spacing, perturb, diagonal swapping). It generates a simple triangulation of a sector, and then use it to generate 6 different triangular grids in a disk or a 50-yen domain. A tecplot style file is included to generate a picture such as the one below after the code is run.

HCH Grid Generation Code Package:   Download a package: hch_v2p1_release.tar.gz This package contains two f90 codes. (1)Structured/unstructured grid generation code for a hemisphere-cylinder-hemisphere body. It generates hex-prism structured, prismatic, tetrahedral, mixed prism-tetra, and mixed prism-hex grids. The shape of the body can be defined by a shape function. Four different functions are implemented in the code: a suboff model, a suboff with a hemisphere front, a sine function, and an ellipsoid. (2)A coarsening program. It generates a series of coarse grids from the grid generated in (1). 10-16-2020: Bugs fixed. 07-24-2019: Added .su2 grid. 12-17-2017: Several bugs fixed.

3D Grid Generation for Hemisphere:   hemisphere_grid_v10.f90 A 3D grid is generated for a hemisphere (with two configurations) and written in the UGRID format. A tetrahedarl, prismatic, or mixed grid can be generated with a specified y-plus value. Learn how a 3D unstructured grid is generated over a hemisphere.

3D Bump Grid Generation:   edu3d_bump_v3.f90 input.nml

Educationally-Designed Unstructured 3D (EDU3D) 3D Tetrahedral Grid Generation:   edu3d_tetgrid_cube_v4.f90 (regular grid) edu3d_tetgrid_cube_ptb_v7.f90 (irregular grid with nodal perturbation) Tetrahedral grid in a unit cube is generated and written in the UGRID format. Learn how the hexahedron can be divided into six tetrahedra. Mdify the code to divide the hexahedron into five tetrahedra. 07-24-19: Added .su2 grid file. 11-24-15: Added an irregular grid version (tetgrid_cube_ptb_v5.f90). 10-06-17: Added a subroutine that computes edge-based data used in the edge-based discretization (tetgrid_cube_ptb_v6p1.f90).

3D Mixed Grid Generation:   mixgrid_cube_v5.f90 Mixed (tetrahedral-prismatic) grid in a unit cube is generated and written in the UGRID format. Learn how a typical viscous-type grid can be generated. You may want to modify the code to apply stretching to the prismatic layer for a smooth transition to the isotropic tetrahedral region. (08-22-2018: A bug fixed about the boundary grid information.)

3D Prismatic Grid Generation:   przgrid_cube_v4.f90 Prismatic grid in a unit cube is generated and written in the UGRID format. It's just a cubic domain, but may be useful in learning unstructured grid data.

3D Hex Grid Generation (Unstructured Format):   hexgrid_cube_v4.f90 Hexahedral grid in a unit cube is generated and written as unstructured/finite-element data in a UGRID 3D unstructured grid file. It may be useful for those who want to learn a typical data structure of unstructured grids.

EDU1D Sock Tube Code:   edu1d_euler_shock_tube.zip This is a code that solves shock tube problems with various flux functions and limiters. It can be useful for learning flux functions, limiters, Hancock's preedictor-corrector scheme, primitive vesus characterisstic limiting, etc., and also useful for computing an exact solution to a shock tube problem. Uploaded. (05-19-19)

How to compute flux Jacobians:   edu3d_exact_flux_jacobian_v0.f90   module_ddt5.f90 [To compile: gfortran edu3d_exact_flux_jacobian_v0.f90] This code shows two ways to compute flux Jacobians: automatic differentiation and hand-coded subroutines for the Roe flux (inviscid) and the alpha-damping flux (viscous) in three dimensions. This serves also as verification of the hand-coded Jacobians (which are complicated). See how the Roe flux is exactly differentiated. It may be useful also for learning how the viscous flux is computed. Uploaded. (04-07-19)

Machine zero:   machinezero_v3.f90 What is machine zero? Maybe, it is a non-zero number which cannot be recognized by a machine. So, it depends on machines. Download, compile, and run it to find out the `zero' in your machine. Comments corrected. (09-23-18) Generalized. Now it finds `machine zero' for a given value. (11-29-12)

Educationally-Designed Unstructured 1D (EDU1D) First/Second/Third-Order Diffusion Schemes in 1D:   edu1d_oned_first_order_diffusion_v1.f90 edu1d_oned_second_order_diffusion_v1.f90 edu1d_oned_third_order_diffusion_v1.f90 Diffusion equation is solved by 1st/2nd/3rd-order upwind schemes on irregularly-spaced grids. The idea is to integrate an equivalent hyperbolic system toward a steady state. This way, we can advance in pseudo time with a large O(h) time step (not O(h^2)), and compute the solution gradient with the equal order of accuracy on irregular grids. Compare with a conventional diffusion scheme to see how accurate the gradient can be and how fast it converges. Note: This third-order scheme does not require computation/storage of second derivatives. See Ref.[2]. References: [1] JCP2014 | Preprint - 1st/2nd/3rd-order FV upwind schemes for diffusion [2] JCP2017 | Preprint - Third-order accuracy without second derivatives.

Solver Package: GCR solver module:     gcr_module_v3.f90 This module is a GCR (Generalized Conjugate Residual) linear solver, which can be used to construct a Jacobian-Free Newton-Krylov solver (JFNK). The code illustrates how a JFNK solver can be implemented.

Educationally-Designed Unstructured 2D (EDU2D) EDU2D-AdvDiff (16 files)   edu2d_advdiff_basic_package_v0.f90   edu2d_advdiff_solver_nc_v2.f90 <--- Corrected   edu2d_advdiff_main_v1.f90   edu2d_advdiff_lsq_grad_nc_v0.f90   edu2d_ddt3_v0.f90   generate_grids_for_edu2d_steady_v0.f90   generate_grids_for_edu2d_unsteady_v0.f90   Makefile   readme_v0.txt   steady_screen_out.txt   unsteady_screen_out.txt   steady_plot_solution.lay   steady_plot_boundary_solution.lay   unsteady_initial_solution.lay   unsteady_final_solution.lay   unsteady_initial_final_solutions.lay Download all files: edu2d-advdiff.zip This is an unstructured advection-diffusion solver. - Node-centered finite-volume discretization - Fully unstructured grids (triangles, quads, or mixed) - Least-squares gradients (linear or quadratic LSQ) - Diffusion is discretized as a hyperbolic system for accurate gradient prediction - Implicit time-stepping with BDF2 - Jacobian-Free Newton-Krylov with GCR - Defect-correction implicit solver as a variable preconditioner - Designed to solve both steady and unsteady problems. Unseady advection diffusion problem. Highly accurate normal gradient prediction on a boundary (steady).

Kahan's sum to minimize the round-off error Avoiding round-off errors:   kahan_sum_example.f90 This code shows how the round-off error can be minimized in computing a sum of small values, by using Kahan's round-off-corrected summation.

Educationally-Designed Unstructured 2D (EDU2D) EDU2D_MMS_TWODNS   edu2d_mms_twodns_v2.f90 This code shows how to compute the source terms in the method of manufactured solutions (MMS) for the 2D Navier-Stokes equations. Any function can be made an exact solution to the 2D Navier-Stokes equations with suitable source terms. That is, any function v(x,y) is an exact solution to the following equation: MMS has been widely used for verifying the order of accuracy of a CFD code. Grab this code, and learn how it can be implemented. I always use MMS to verify the implementation. It is a very useful tool for debugging.

Educationally-Designed Unstructured 2D (EDU2D) EDU2D_Quadrature_Triangle   edu2d_quadrature_triangle_v1.f90 This program illustrates the use of a Gaussian quadrature, which integrates a function over a triangle. An example function, x^5 + y^5 -32*x^2*y^3 + x^3*y^2, is integrated over a triangle defined by the vertices: (0.1,0.3), (1.3,-0.2), (-0.9,1.7), and the numerical value obtained by a 7-point quadrature formula is compared with the exact value 0.53726127666666666668, up to 15 significant digits (double precision).

Educationally-Designed Unstructured 2D (EDU2D) EDU2D-Airfoil-Spline (3 files)   edu2d_airfoil_spline_v1.f90   om6_wing_section_sharp.dat   edu2d_airfoil_plot.m Compile edu2d_airfoil_spline_v1.f90 and run it. This f90 program reads OM6 wing section data contained in the file 'om6_wing_section_sharp.dat', constructs a cubic spline (piecewise cubic polynomial interpolation), and writes out data files for plotting the airfoil by Matlab or Tecplot. I remember I used a cubic spline to represent a Joukowsky airfoil, so that boundary nodes can be adaptively moved along the surface to better represent numerical solutions [ IJNMF2002 ] The OM6 section data are from [ AIAA J. Vol. 54, No. 9, September 2016 ], and the trailing edge is closed as described in the AIAA jounal paper.

Educationally-Designed Unstructured 2D (EDU2D) EDU2D-Tecplot2Grid (5 files)   edu2d_tecplot2grid_v1.f90   sample_input   tria_rankine_grid_tecplot.dat   tria_rankine_grid_tecplot.data   readme_v1.txt Download all files: edu2d_tecplot2grid.zip This program reads a Tecplot data file, generates a boundary information, and writes .grid and .bcmap file for edu2d solvers. Tecplot file does not contain boundary information, and so the boundary information is genearted inside the program. Not an elegant code, but it may be useful to learn how to identify boundary nodes for a given grid, which is not really trivial. A grid for a Rankine's half body is used as an example.

Educationally-Designed Unstructured 2D (EDU2D) EDU2D-Euler-Steady (9 files)   edu2d_euler_steady_basic_package_v0.f90   edu2d_euler_steady_solver_v0.f90   edu2d_euler_steady_main_v0.f90   edu2d_euler_jacobian_v1.f90   edu2d_euler_linear_solve_v0.f90   twod_bump_irregular_grid_v0.f90   Makefile   readme_v0.txt   bump_screen_out.txt Download all files: edu2d-euler-steady.zip Want to learn how to write an implicit unstructured CFD code? Grab this code, look inside to see how it is written, get good understanding, and then write your own. This code computes a steady flow over a bump with the Roe flux by two solution methods: an explicit 2-stage Runge-Kutta scheme and an implicit (defect correction) method with the exact Jacobian for a 1st-order scheme, on irregular triangular grids. A grid generation code is included for a bump problem. - Node-centered finite-volume discretization - Fully unstructured grids (triangles, quads, or mixed) - Least-squares gradients (linear or quadratic LSQ) - Defect-correction implicit solver (Newton's method for 1st-order) - Both explicit and implicit solvers are implemented

Educationally-Designed Unstructured 2D (EDU2D) EDU2D-Euler-Unsteady (8 files)   edu2d_euler_unsteady_basic_package_v0.f90   edu2d_euler_unsteady_solver_v0.f90   edu2d_euler_unsteady_main_v0.f90   twod_rectangular_grid_v0.f90   Makefile   readme_v0.txt   project_quad_screen_out.txt   project_tria_screen_out.txt Download all files: edu2d-euler-unsteady.zip Want to learn how to write an unstructured CFD code? Grab this code, look inside to see how it is written, get good understanding, and then write your own. This code has Roe and Rotated-RHLL fluxes, Van Albada limiter, and a 2-stage Runge-Kutta time-stepping for solving a shock diffraction problem. It works for quadrilateral grids, triangular grids, and mixed grids also. It is set up to solve a shock diffraction problem. You can easily modify it for solving other problems. A grid generation code is included. - Node-centered finite-volume discretization - Fully unstructured grids (triangles, quads, or mixed) - Least-squares gradients (linear or quadratic LSQ) - Limiter is applied at edges - Explicit time-stepping with 2-stage Runge-Kutta scheme

Educationally-Designed Unstructured 2D (EDU2D) EDU2D-Template: 2D Example Unstructured CFD Code (8 files)   edu2d_basic_package_v0.f90   edu2d_template_main_v0.f90   edu2d_template_nc_v0.f90   edu2d_template_cc_v0.f90   generate_grids_for_edu2d_v0.f90   Makefile   readme_v0.txt   any_name_screen_out.txt Download all files: edu2d-template.zip This code reads an unstructured grid file, generate various grid data, go through some dummy CFD solvers, and then writes out Tecplot data files for viewing the solution and the grid. Solvers are dummies and so do not solve anything, but you'll see how a node/cell-centered finite-volume solver can be implemented, e.g., loop over edges, nodes, elements, and computation of least-squares gradients, etc. You can generate your own unstructured CFD code by replacing the dummy solver by your own solver. A grid generation file is included, which can be used to generate grid files you need to run the EDU2D-Template code.

MPI program for Jacobi Iteration Ver.1:   jacobi_mpi_v1.f90 This is an MPI program for the Jacobi iteration, solving the finite-difference discretization of the Laplace equation in a square domain. It can be useful for learning how to write an MPI program.

Automatic Differentiation Ver.1:   ad_driver_v1.f90   module_ddt3.f90 Automatic differentiation is nice. It is sometimes used to compute flux Jacobians for implicit formulaitons in CFD. To those who are curious about it, here is a set of files with which you can experience and learn how it works. Simply compile the file "ad_driver.f90" and run it (the other file will be automatically included in the program) to get a feeling of automatic differentiation.

Numerical Fluxes (3D Euler) Ver.3: threed_euler_fluxes_v3.f90 Here are 3D Euler numerical fluxes. Download and take a look. Learn how the standard Roe and a very robust rotated-hybrid fluxes are implemented for the 3D Euler equations, and also how the Roe flux can be implemented without tangent vectors. Included are Roe with/without tangent vectors, and Rotated-RHLL fluxes. I'll be extremely happy if you kindly report bugs. Thank you, arigatou! Bugs fixed for the RHLL flux, more comments added, Roe without tangent vectors added. (04-30-12).

3D Grid Generation for Hemishpere Cylinder: The code is available at NASA's Turbulence Modeling Resource (TMR) website.

1D Hyperbolic Diffusion Scheme:   oned_upwind_diffusion.f90 Believe or not, the diffusion equation is solved by an upwind scheme. The idea is to integrate an equivalent hyperbolic system toward a steady state. This way, we can advance in time with a large O(h) time step (not O(h^2)), and compute the solution gradient with the equal order of accuracy. Compare with a common scheme (Galerkin) for 512 nodes to see how fast the upwind scheme can be. Reference: JCP2007 | Preprint

Genetic Algorithm, Ver.1:   ga_v1.f90 This is a simple genetic algorithm program for finding minimizers of a function. I wrote this in 1999 when I was interested in aerodynamic optimization problems. I think that genetic algorithm is a very interesting optimization algorithm.

1D Euler Code Ver.1:   oned_euler_v1.f90 ,   oned_euler_plot_v1.m Here is a 1D Euler code (1D shock tube code) for solving Sod's shock tube problem, using Roe's Approximate Riemann solver, minmod limiter, and 2-stage Runge-Kutta time-stepping. Learn how a second-order non-oscillatory Euler code is written, or just run it to see how it is capable of computing discontinuous solutions. Incorporate various flux subroutines given below to explore other methods. Minor bugs fixed (12-29-10).

Numerical Fluxes (1D Euler) Ver.5:   oned_euler_fluxes_v5.f90 All numerical fluxes are functions of the left and right states (and possibly dt and dx). Learn how those famous fluxes can be implemented, or just use them to see how they work for various shock-tube problems. Included are Lax-Friedrichs, Richtmyer, MacCormack, Steger-Warming, Van Leer, AUSM, Zha-Bilgen, Godunov, Osher, Roe, Rusanov, HLL, HLLL, AUFS flux functions. Bug fixed for Godunov fluxes.

Numerical Fluxes (2D Euler) Ver.2:   twod_euler_fluxes_v2.f90 Want 2D fluxes also? Here you are. Download and take a look. Learn how the standard Roe and a very robust rotated-hybrid fluxes are implemented. Included are Roe and Rotated-RHLL fluxes (Riemann Solvers). Bugs fixed for the Rotated-RHLL flux. It works very well now.

Blasius solutions for a flow over a flat-plate (compute the exact solution) :   blasius_v1.f90 ,   blasius_plot_v1.m Exact solution for a flow over a flat-plate. This solution has been used by many people to verify the accuracy of their Navier-Stokes code. If you have not, use this code now to compute the exact solution at any point on your grid, compute the error, and verify the accuracy of your code. Minor bugs fixed (12-29-10).

Viscous shock structure NS-solution (compute the exact solution) :   ns_shock_structure_v1.f90 ,   ns_shock_structure_plot_v1.m Exact solution for a shock wave internal structure to the 1D Navier-Stokes equations. This solution has been used by some people to verify the accuracy of their 1D Navier-Stokes code. If you have not, use this code now to generate a data file for the exact solution, use it as an initial solution for your code, converge to a steady state, and verify the accuracy of your code. Minor bugs fixed (05-17-2023, 12-29-10).