diff --git a/.appveyor.yml b/.appveyor.yml
index 11a0635..f091cc5 100644
--- a/.appveyor.yml
+++ b/.appveyor.yml
@@ -9,6 +9,7 @@ platform: x64
environment:
PYTHON: C:\Python35-x64
PATH: C:\Program Files\Project\bin;C:\Program Files\Microsoft MPI\Bin;%PATH%
+ AITHER_FLUID_DATABASE: C:\Program Files\Project\fluidDatabase
build_script:
- ps: >-
powershell ci\appveyor\installMPI.ps1
@@ -18,4 +19,4 @@ test_script:
- cmd: >-
cd testCases
- %PYTHON%\python.exe regressionTests.py --mpirunPath=mpiexec.exe --aitherPath=aither.exe
\ No newline at end of file
+ %PYTHON%\python.exe regressionTests.py --mpirunPath=mpiexec.exe --aitherPath=aither.exe --operatingSystem=windows
\ No newline at end of file
diff --git a/.github/ISSUE_TEMPLATE.md b/.github/ISSUE_TEMPLATE.md
new file mode 100644
index 0000000..5f019ac
--- /dev/null
+++ b/.github/ISSUE_TEMPLATE.md
@@ -0,0 +1,11 @@
+## Issue Type (bug, feature request, other)
+
+
+## Description of Issue
+
+
+## Steps to Reproduce Issue
+
+
+## System Details (OS, compiler, number of processors, etc)
+
diff --git a/.github/PULL_REQUEST_TEMPLATE.md b/.github/PULL_REQUEST_TEMPLATE.md
new file mode 100644
index 0000000..f2eef47
--- /dev/null
+++ b/.github/PULL_REQUEST_TEMPLATE.md
@@ -0,0 +1,8 @@
+## Fixes Issue #
+
+
+## Description of Changes
+
+
+## Suggested Reviewers
+@mnucci32
diff --git a/.travis.yml b/.travis.yml
index 142c6af..77e4540 100644
--- a/.travis.yml
+++ b/.travis.yml
@@ -27,6 +27,7 @@ matrix:
- CXX_COMPILER=g++-5
- C_COMPILER=gcc-5
- BUILD_TYPE=release
+ - AITHER_FLUID_DATABASE=$TRAVIS_BUILD_DIR/build/install/fluidDatabase
- os: linux
dist: trusty
sudo: required
@@ -41,6 +42,22 @@ matrix:
- CXX_COMPILER=g++-6
- C_COMPILER=gcc-6
- BUILD_TYPE=debug
+ - AITHER_FLUID_DATABASE=$TRAVIS_BUILD_DIR/build/install/fluidDatabase
+ - os: linux
+ dist: trusty
+ sudo: required
+ compiler: gcc
+ addons:
+ apt:
+ sources:
+ - ubuntu-toolchain-r-test
+ packages:
+ - g++-7 gcc-7 libstdc++-7-dev
+ env:
+ - CXX_COMPILER=g++-7
+ - C_COMPILER=gcc-7
+ - BUILD_TYPE=release
+ - AITHER_FLUID_DATABASE=$TRAVIS_BUILD_DIR/build/install/fluidDatabase
- os: linux
dist: trusty
sudo: required
@@ -56,6 +73,7 @@ matrix:
- CXX_COMPILER=clang++-3.8
- C_COMPILER=clang-3.8
- BUILD_TYPE=release
+ - AITHER_FLUID_DATABASE=$TRAVIS_BUILD_DIR/build/install/fluidDatabase
- os: osx
osx_image: xcode8
compiler: gcc
@@ -65,6 +83,7 @@ matrix:
- HOMEBREW_CC=gcc-7
- HOMEBREW_CXX=g++-7
- BUILD_TYPE=release
+ - AITHER_FLUID_DATABASE=$TRAVIS_BUILD_DIR/build/install/fluidDatabase
- os: osx
osx_image: xcode8
compiler: clang
@@ -74,6 +93,7 @@ matrix:
- HOMEBREW_CC=clang
- HOMEBREW_CXX=clang++
- BUILD_TYPE=release
+ - AITHER_FLUID_DATABASE=$TRAVIS_BUILD_DIR/build/install/fluidDatabase
# upgrade packages
before_install:
@@ -100,10 +120,10 @@ before_script:
# build instructions
script:
- cmake -G "Unix Makefiles" -DCMAKE_CXX_COMPILER=$CXX_COMPILER -DCMAKE_C_COMPILER=$C_COMPILER -DMPI_DIR=$TRAVIS_BUILD_DIR/openmpi -DCMAKE_BUILD_TYPE=$BUILD_TYPE -DCMAKE_INSTALL_PREFIX=install ..
- - make -j2
+ - make -j 4
- make install
- cd ../testCases
- - python3 regressionTests.py --aitherPath=$TRAVIS_BUILD_DIR/build/install/bin/aither --operatingSystem=$TRAVIS_OS_NAME --mpirunPath=$TRAVIS_BUILD_DIR/openmpi/bin/mpirun
+ - python3 regressionTests.py --aitherPath=$TRAVIS_BUILD_DIR/build/install/bin/aither --operatingSystem=$TRAVIS_OS_NAME --mpirunPath=$TRAVIS_BUILD_DIR/openmpi/bin/mpirun --build $BUILD_TYPE
- cd $TRAVIS_BUILD_DIR/build
after_success:
diff --git a/CODE_OF_CONDUCT.md b/CODE_OF_CONDUCT.md
new file mode 100644
index 0000000..61814c7
--- /dev/null
+++ b/CODE_OF_CONDUCT.md
@@ -0,0 +1,46 @@
+# Contributor Covenant Code of Conduct
+
+## Our Pledge
+
+In the interest of fostering an open and welcoming environment, we as contributors and maintainers pledge to making participation in our project and our community a harassment-free experience for everyone, regardless of age, body size, disability, ethnicity, gender identity and expression, level of experience, nationality, personal appearance, race, religion, or sexual identity and orientation.
+
+## Our Standards
+
+Examples of behavior that contributes to creating a positive environment include:
+
+* Using welcoming and inclusive language
+* Being respectful of differing viewpoints and experiences
+* Gracefully accepting constructive criticism
+* Focusing on what is best for the community
+* Showing empathy towards other community members
+
+Examples of unacceptable behavior by participants include:
+
+* The use of sexualized language or imagery and unwelcome sexual attention or advances
+* Trolling, insulting/derogatory comments, and personal or political attacks
+* Public or private harassment
+* Publishing others' private information, such as a physical or electronic address, without explicit permission
+* Other conduct which could reasonably be considered inappropriate in a professional setting
+
+## Our Responsibilities
+
+Project maintainers are responsible for clarifying the standards of acceptable behavior and are expected to take appropriate and fair corrective action in response to any instances of unacceptable behavior.
+
+Project maintainers have the right and responsibility to remove, edit, or reject comments, commits, code, wiki edits, issues, and other contributions that are not aligned to this Code of Conduct, or to ban temporarily or permanently any contributor for other behaviors that they deem inappropriate, threatening, offensive, or harmful.
+
+## Scope
+
+This Code of Conduct applies both within project spaces and in public spaces when an individual is representing the project or its community. Examples of representing a project or community include using an official project e-mail address, posting via an official social media account, or acting as an appointed representative at an online or offline event. Representation of a project may be further defined and clarified by project maintainers.
+
+## Enforcement
+
+Instances of abusive, harassing, or otherwise unacceptable behavior may be reported by contacting the project team at michael.nucci@gmail.com. The project team will review and investigate all complaints, and will respond in a way that it deems appropriate to the circumstances. The project team is obligated to maintain confidentiality with regard to the reporter of an incident. Further details of specific enforcement policies may be posted separately.
+
+Project maintainers who do not follow or enforce the Code of Conduct in good faith may face temporary or permanent repercussions as determined by other members of the project's leadership.
+
+## Attribution
+
+This Code of Conduct is adapted from the [Contributor Covenant][homepage], version 1.4, available at [http://contributor-covenant.org/version/1/4][version]
+
+[homepage]: http://contributor-covenant.org
+[version]: http://contributor-covenant.org/version/1/4/
diff --git a/CONTRIBUTING.md b/CONTRIBUTING.md
new file mode 100644
index 0000000..910a713
--- /dev/null
+++ b/CONTRIBUTING.md
@@ -0,0 +1,40 @@
+# How to Contribute
+Thank you for contributing to aither! Check out the aither [blog](http://aithercfd.com) for more information. Aither
+uses the git flow branching method. The **master** branch is a permanent branch containing stable code with releases
+coming off of this branch. The default branch is the **develop** branch. This is a permanent branch containing the most
+up-to-date code and is therefore not guaranteed to be working. Bug fixes and feature additions branch off of **develop** and
+are merged back in when complete. Aither uses continuous integration services [TravisCI](https://travis-ci.org/mnucci32/aither)
+and [Appveyor](https://ci.appveyor.com/project/mnucci32/aither/branch/develop) to test builds on Linux, macOS, and
+Windows. The code is built on each platform and a suite of tests are run to ensure that the residuals have not
+changed due to a commit. There are a few ways to get involved with this project including reporting or fixing a bug,
+requesting or adding a feature, and adding documentation.
+
+### Reporting a Bug
+Bugs should be reported via [Github Issues](https://github.com/mnucci32/aither/issues).
+
+### Fixing a Bug
+A list of known bugs can be found on the [Issues](https://github.com/mnucci32/aither/issues) page. A fix can be submitted
+by creating a new branch off of the **develop** branch with a descriptive name starting with **hotfix_**. When the bug has
+been fixed and and all continuous integration tests are passing, a pull request can be submitted to merge the fix back into
+**develop**.
+
+### Requesting a Feature
+New features should be requested via [Github Issues](https://github.com/mnucci32/aither/issues). Not all requests will be
+accepted.
+
+### Adding a Feature
+A list of features to be added can be found on the [Issues](https://github.com/mnucci32/aither/issues) page. If you want
+to add a feature not listed, first submit an issue requesting that it be added. New features should be implemented by
+creating a branch off of **develop** with a descriptive name starting with **feature_**. When the feature has been fixed
+and all continuous integreation tests are passing, a pull request can be submitted to merge the feature back into **develop**.
+If it is a fairly significant addition, an additional test case should be added as well.
+
+### Adding Documentation
+Documentation can be added by editing the [Wiki](https://github.com/mnucci32/aither/wiki) page.
+
+# Coding Style
+* Aither conforms to the [Google Style Guide](https://google.github.io/styleguide/cppguide.html)
+* Variable names should use camelCase starting with a lower case letter
+* Function names should use CamelCase starting with an upper case letter
+* Class member variables should be appended with a _
+* Modern C++ features should be used when possible
diff --git a/Makefile b/Makefile
deleted file mode 100644
index ac8bf0d..0000000
--- a/Makefile
+++ /dev/null
@@ -1,77 +0,0 @@
-OBJS = main.o plot3d.o input.o boundaryConditions.o eos.o primVars.o procBlock.o output.o matrix.o parallel.o slices.o turbulence.o inviscidFlux.o viscousFlux.o source.o resid.o kdtree.o genArray.o fluxJacobian.o uncoupledScalar.o utility.o
-CC = mpic++
-DEBUG = -O0 -ggdb -pg
-OPTIM = -O3 -march=native
-PROF = -O3 -march=native -pg
-CODENAME = aither
-CFLAGS = -std=c++14 -Wall -pedantic -c $(OPTIM)
-LFLAGS = -std=c++14 -Wall -pedantic $(OPTIM) -o $(CODENAME)
-
-$(CODENAME) : $(OBJS)
- $(CC) $(LFLAGS) $(OBJS)
-
-plot3d.o : plot3d.cpp plot3d.hpp vector3d.hpp multiArray3d.hpp
- $(CC) $(CFLAGS) plot3d.cpp
-
-main.o : main.cpp plot3d.hpp vector3d.hpp input.hpp procBlock.hpp eos.hpp primVars.hpp boundaryConditions.hpp inviscidFlux.hpp tensor.hpp viscousFlux.hpp output.hpp parallel.hpp turbulence.hpp resid.hpp multiArray3d.hpp genArray.hpp fluxJacobian.hpp utility.hpp
- $(CC) $(CFLAGS) main.cpp
-
-input.o : input.cpp input.hpp boundaryConditions.hpp
- $(CC) $(CFLAGS) input.cpp
-
-primVars.o : primVars.cpp primVars.hpp vector3d.hpp eos.hpp inviscidFlux.hpp boundaryConditions.hpp input.hpp macros.hpp genArray.hpp
- $(CC) $(CFLAGS) primVars.cpp
-
-procBlock.o : procBlock.cpp procBlock.hpp vector3d.hpp plot3d.hpp eos.hpp primVars.hpp inviscidFlux.hpp input.hpp genArray.hpp viscousFlux.hpp boundaryConditions.hpp macros.hpp turbulence.hpp kdtree.hpp uncoupledScalar.hpp fluxJacobian.hpp matrix.hpp utility.hpp
- $(CC) $(CFLAGS) procBlock.cpp
-
-inviscidFlux.o : inviscidFlux.cpp vector3d.hpp eos.hpp primVars.hpp inviscidFlux.hpp input.hpp macros.hpp genArray.hpp turbulence.hpp matrix.hpp
- $(CC) $(CFLAGS) inviscidFlux.cpp
-
-boundaryConditions.o : boundaryConditions.cpp boundaryConditions.hpp plot3d.hpp vector3d.hpp
- $(CC) $(CFLAGS) boundaryConditions.cpp
-
-eos.o : eos.cpp eos.hpp vector3d.hpp
- $(CC) $(CFLAGS) eos.cpp
-
-slices.o : slices.cpp slices.hpp procBlock.hpp
- $(CC) $(CFLAGS) slices.cpp
-
-viscousFlux.o : viscousFlux.cpp vector3d.hpp tensor.hpp eos.hpp primVars.hpp viscousFlux.hpp input.hpp turbulence.hpp macros.hpp
- $(CC) $(CFLAGS) viscousFlux.cpp
-
-output.o : output.cpp output.hpp procBlock.hpp tensor.hpp vector3d.hpp plot3d.hpp eos.hpp primVars.hpp inviscidFlux.hpp input.hpp turbulence.hpp genArray.hpp
- $(CC) $(CFLAGS) output.cpp
-
-parallel.o : parallel.cpp parallel.hpp primVars.hpp procBlock.hpp vector3d.hpp plot3d.hpp boundaryConditions.hpp resid.hpp
- $(CC) $(CFLAGS) parallel.cpp
-
-matrix.o : matrix.cpp matrix.hpp macros.hpp genArray.hpp
- $(CC) $(CFLAGS) matrix.cpp
-
-genArray.o : genArray.cpp genArray.hpp macros.hpp
- $(CC) $(CFLAGS) genArray.cpp
-
-turbulence.o : turbulence.cpp turbulence.hpp matrix.hpp vector3d.hpp tensor.hpp primVars.hpp eos.hpp
- $(CC) $(CFLAGS) turbulence.cpp
-
-source.o : source.cpp source.hpp macros.hpp turbulence.hpp primVars.hpp matrix.hpp
- $(CC) $(CFLAGS) source.cpp
-
-resid.o : resid.cpp resid.hpp
- $(CC) $(CFLAGS) resid.cpp
-
-kdtree.o : kdtree.cpp kdtree.hpp vector3d.hpp
- $(CC) $(CFLAGS) kdtree.cpp
-
-fluxJacobian.o : fluxJacobian.cpp fluxJacobian.hpp turbulence.hpp vector3d.hpp primVars.hpp eos.hpp input.hpp genArray.hpp matrix.hpp inviscidFlux.hpp uncoupledScalar.hpp tensor.hpp utility.hpp
- $(CC) $(CFLAGS) fluxJacobian.cpp
-
-uncoupledScalar.o : uncoupledScalar.cpp uncoupledScalar.hpp genArray.hpp
- $(CC) $(CFLAGS) uncoupledScalar.cpp
-
-utility.o : utility.cpp utility.hpp genArray.hpp vector3d.hpp multiArray3d.hpp procBlock.hpp eos.hpp input.hpp turbulence.hpp slices.hpp fluxJacobian.hpp kdtree.hpp resid.hpp
- $(CC) $(CFLAGS) utility.cpp
-
-clean:
- rm *.o *~ $(CODENAME)
diff --git a/README.md b/README.md
index e23163e..8093bca 100644
--- a/README.md
+++ b/README.md
@@ -5,11 +5,11 @@
| Master | [![Build Status](https://travis-ci.org/mnucci32/aither.svg?branch=master)](https://travis-ci.org/mnucci32/aither) | [![Build status](https://ci.appveyor.com/api/projects/status/o7fc231lp9jxlsib/branch/master?svg=true)](https://ci.appveyor.com/project/mnucci32/aither/branch/master) | [![Coverage Status](https://codecov.io/github/mnucci32/aither/coverage.svg?branch=master)](https://codecov.io/github/mnucci32/aither?branch=master) |
| Develop | [![Build Status](https://travis-ci.org/mnucci32/aither.svg?branch=develop)](https://travis-ci.org/mnucci32/aither) | [![Build status](https://ci.appveyor.com/api/projects/status/o7fc231lp9jxlsib/branch/develop?svg=true)](https://ci.appveyor.com/project/mnucci32/aither/branch/develop) | [![Coverage Status](https://codecov.io/github/mnucci32/aither/coverage.svg?branch=develop)](https://codecov.io/github/mnucci32/aither?branch=develop) |
-### About The code
+### About The Code
This code is for a 3D Navier-Stokes computational fluid dynamics solver. It is
a cell centered, structured solver, using multi-block structured grids in Plot3D
format. It uses explicit and implicit time integration methods. It uses MUSCL
-extrapolation to reconstruct the primative variables from the cell centers to
+extrapolation to reconstruct the primitive variables from the cell centers to
the cell faces for 2nd order accuracy. Higher order reconstruction is acheived
with a 5th order WENO reconstruction for the inviscid fluxes, and a 4th order
central reconstruction for the viscous fluxes. The code uses the Roe
@@ -28,31 +28,39 @@ are the implicit euler (1st order), Crank-Nicholson (2nd order), and BDF2
parallel using MPI. For RANS simulations the Wilcox K-Omega 2006 and SST 2003
turbulence models are available. Wall functions are supported for both models.
For detatched eddy simulations, the SST-DES turbulence model is available. For
-large eddy simulations, the WALE subgrid scale model is available.
+large eddy simulations, the WALE subgrid scale model is available. A
+Multi-species flow capability is in progress. Supported diffusion models will be
+Schmidt number based diffusion.
### To Do List
-* Add non-reflecting boundary conditions
* Add multigrid scheme for improved convergence
-* Add multi-species flow capability
+* Add reacting flow capability
+* Performance improvements
### Dependencies
* MPI - OpenMPI, MPICH, & MS-MPI have been used
* C++ compiler with C++14 support
* Cmake - Cmake only depends on a C++ compiler
-### How To compile
+### How To Compile And Install
Aither is compiled and installed with the standard cmake process.
```bash
cmake -DCMAKE_INSTALL_PREFIX=/path/to/installation -DCMAKE_BUILD_TYPE=release /path/to/source
make
make install
+
+export AITHER_FLUID_DATABASE=/path/to/installation/fluidDatabase
```
Cmake will automatically look for an MPI package. To specify a specific
-installation, set *-DMPI_DIR* to the MPI installation directory. In addition
-to *release*, other supported build types are *debug*, *profile*,
-*relwithdebinfo*, and *minsizerel*.
+installation, set **-DMPI_DIR** to the MPI installation directory. In addition
+to **release**, other supported build types are **debug**, **profile**,
+**relwithdebinfo**, and **minsizerel**.
+
+The **AITHER_FLUID_DATABASE** environment variable should be set to the
+**fluidDatabase** folder inside of the installation directory. This tells aither
+where to look for fluid properties.
### How To Run
```bash
diff --git a/ci/appveyor/buildAither.ps1 b/ci/appveyor/buildAither.ps1
index e221507..1da5ca5 100644
--- a/ci/appveyor/buildAither.ps1
+++ b/ci/appveyor/buildAither.ps1
@@ -4,7 +4,7 @@ function BuildAither() {
md build
cd build
cmake -G "Visual Studio 14 2015 Win64" -DMPI_DIR="C:\Program Files (x86)\Microsoft SDKs\MPI" ..
- cmake --build . --target INSTALL --config Release
+ cmake --build . --target INSTALL --config release
cd ..
}
diff --git a/ci/appveyor/installMPI.ps1 b/ci/appveyor/installMPI.ps1
index d2046db..5dfc293 100644
--- a/ci/appveyor/installMPI.ps1
+++ b/ci/appveyor/installMPI.ps1
@@ -9,11 +9,24 @@ function InstallMPI() {
Write-Host "Installing Microsoft MPI Runtime..."
appveyor DownloadFile http://download.microsoft.com/download/B/2/E/B2EB83FE-98C2-4156-834A-E1711E6884FB/MSMpiSetup.exe
Start-Process -FilePath MSMpiSetup.exe -ArgumentList -unattend -Wait
+ Write-Host "Microsoft MPI Runtime installation complete..."
cd ..
}
+function InstallCmake() {
+ Write-Host "Installing CMake 3.4.0 ..." -ForegroundColor Cyan
+ $exePath = "$($env:USERPROFILE)\cmake-3.4.0-rc2-win32-x86.exe"
+ Write-Host "Downloading..."
+ (New-Object Net.WebClient).DownloadFile('https://cmake.org/files/v3.4/cmake-3.4.0-rc2-win32-x86.exe', $exePath)
+ Write-Host "Installing..."
+ cmd /c start /wait $exePath /S
+ cmake --version
+ Write-Host "CMake 3.4.0 installed" -ForegroundColor Green
+}
+
function main() {
InstallMPI
+ InstallCmake
}
main
\ No newline at end of file
diff --git a/fluidDatabase/Ar.dat b/fluidDatabase/Ar.dat
new file mode 100644
index 0000000..a3f1412
--- /dev/null
+++ b/fluidDatabase/Ar.dat
@@ -0,0 +1,17 @@
+# general fluid properties
+n: 1.5
+molarMass: 39.948
+vibrationalTemperature: []
+
+# Sutherland's coefficients for transport properties
+# ----- viscosity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandViscosityC1: 2.0343E-06
+sutherlandViscosityS: 1.6053E+02
+
+# ----- thermal conductivity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandConductivityC1: 1.5877E-03
+sutherlandConductivityS: 1.6053E+02
\ No newline at end of file
diff --git a/fluidDatabase/CH4.dat b/fluidDatabase/CH4.dat
new file mode 100644
index 0000000..623df33
--- /dev/null
+++ b/fluidDatabase/CH4.dat
@@ -0,0 +1,17 @@
+# general fluid properties
+n: 3.0
+molarMass: 16.0425
+vibrationalTemperature: [4196.38, 2207.18, 2207.18, 4343.43, 4343.43, 4343.43, 1879.13, 1879.13, 1879.13]
+
+# Sutherland's coefficients for transport properties
+# ----- viscosity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandViscosityC1: 1.0166E-06
+sutherlandViscosityS: 1.6471E+02
+
+# ----- thermal conductivity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandConductivityC1: 1.7680E-02
+sutherlandConductivityS: 2.3083E+03
\ No newline at end of file
diff --git a/fluidDatabase/CO.dat b/fluidDatabase/CO.dat
new file mode 100644
index 0000000..7ef4cfa
--- /dev/null
+++ b/fluidDatabase/CO.dat
@@ -0,0 +1,17 @@
+# general fluid properties
+n: 2.5
+molarMass: 28.0101
+vibrationalTemperature: [3121.5]
+
+# Sutherland's coefficients for transport properties
+# ----- viscosity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandViscosityC1: 1.4500E-06
+sutherlandViscosityS: 1.2882E+02
+
+# ----- thermal conductivity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandConductivityC1: 2.6880E-03
+sutherlandConductivityS: 2.7617E+02
\ No newline at end of file
diff --git a/fluidDatabase/CO2.dat b/fluidDatabase/CO2.dat
new file mode 100644
index 0000000..3ed321a
--- /dev/null
+++ b/fluidDatabase/CO2.dat
@@ -0,0 +1,17 @@
+# general fluid properties
+n: 2.5
+molarMass: 44.0095
+vibrationalTemperature: [960.1, 960.1, 1932.1, 3380.1]
+
+# Sutherland's coefficients for transport properties
+# ----- viscosity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandViscosityC1: 1.6491E-06
+sutherlandViscosityS: 2.6968E+02
+
+# ----- thermal conductivity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandConductivityC1: 4.1247E-03
+sutherlandConductivityS: 8.8020E+02
\ No newline at end of file
diff --git a/fluidDatabase/H.dat b/fluidDatabase/H.dat
new file mode 100644
index 0000000..67e72c5
--- /dev/null
+++ b/fluidDatabase/H.dat
@@ -0,0 +1,17 @@
+# general fluid properties
+n: 1.5
+molarMass: 1.00794
+vibrationalTemperature: []
+
+# Sutherland's coefficients for transport properties
+# ----- viscosity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandViscosityC1: 8.4958E-07
+sutherlandViscosityS: 1.6775E+02
+
+# ----- thermal conductivity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandConductivityC1: 2.6278E-02
+sutherlandConductivityS: 1.6775E+02
\ No newline at end of file
diff --git a/fluidDatabase/H2.dat b/fluidDatabase/H2.dat
new file mode 100644
index 0000000..b7152a7
--- /dev/null
+++ b/fluidDatabase/H2.dat
@@ -0,0 +1,17 @@
+# general fluid properties
+n: 2.5
+molarMass: 2.01588
+vibrationalTemperature: [6338.3]
+
+# Sutherland's coefficients for transport properties
+# ----- viscosity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandViscosityC1: 6.8021E-07
+sutherlandViscosityS: 1.0031E+02
+
+# ----- thermal conductivity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandConductivityC1: 1.5056E-02
+sutherlandConductivityS: 1.3207E+02
\ No newline at end of file
diff --git a/fluidDatabase/H2O.dat b/fluidDatabase/H2O.dat
new file mode 100644
index 0000000..00d09a2
--- /dev/null
+++ b/fluidDatabase/H2O.dat
@@ -0,0 +1,17 @@
+# general fluid properties
+n: 3.0
+molarMass: 18.0153
+vibrationalTemperature: [2294.3, 5261.7, 5403.8]
+
+# Sutherland's coefficients for transport properties
+# ----- viscosity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandViscosityC1: 1.9293E-06
+sutherlandViscosityS: 7.0274E+02
+
+# ----- thermal conductivity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandConductivityC1: 1.1200E-02
+sutherlandConductivityS: 2.0728E+03
\ No newline at end of file
diff --git a/fluidDatabase/He.dat b/fluidDatabase/He.dat
new file mode 100644
index 0000000..2515fde
--- /dev/null
+++ b/fluidDatabase/He.dat
@@ -0,0 +1,17 @@
+# general fluid properties
+n: 1.5
+molarMass: 4.002602
+vibrationalTemperature: []
+
+# Sutherland's coefficients for transport properties
+# ----- viscosity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandViscosityC1: 1.4872E-06
+sutherlandViscosityS: 9.7629E+01
+
+# ----- thermal conductivity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandConductivityC1: 1.1584E-02
+sutherlandConductivityS: 9.7629E+01
\ No newline at end of file
diff --git a/fluidDatabase/N.dat b/fluidDatabase/N.dat
new file mode 100644
index 0000000..3ea3374
--- /dev/null
+++ b/fluidDatabase/N.dat
@@ -0,0 +1,17 @@
+# general fluid properties
+n: 1.5
+molarMass: 14.0067
+vibrationalTemperature: []
+
+# Sutherland's coefficients for transport properties
+# ----- viscosity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandViscosityC1: 1.2953E-06
+sutherlandViscosityS: 1.1190E+02
+
+# ----- thermal conductivity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandConductivityC1: 2.8831E-03
+sutherlandConductivityS: 1.1190E+02
\ No newline at end of file
diff --git a/fluidDatabase/N2.dat b/fluidDatabase/N2.dat
new file mode 100644
index 0000000..bd40d1c
--- /dev/null
+++ b/fluidDatabase/N2.dat
@@ -0,0 +1,17 @@
+# general fluid properties
+n: 2.5
+molarMass: 28.0134
+vibrationalTemperature: [3392]
+
+# Sutherland's coefficients for transport properties
+# ----- viscosity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandViscosityC1: 1.4742E-06
+sutherlandViscosityS: 1.2846E+02
+
+# ----- thermal conductivity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandConductivityC1: 2.6834E-03
+sutherlandConductivityS: 2.5615E+02
\ No newline at end of file
diff --git a/fluidDatabase/NO.dat b/fluidDatabase/NO.dat
new file mode 100644
index 0000000..86c3022
--- /dev/null
+++ b/fluidDatabase/NO.dat
@@ -0,0 +1,17 @@
+# general fluid properties
+n: 2.5
+molarMass: 30.0061
+vibrationalTemperature: [2739]
+
+# Sutherland's coefficients for transport properties
+# ----- viscosity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandViscosityC1: 1.5257E-06
+sutherlandViscosityS: 1.2846E+02
+
+# ----- thermal conductivity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandConductivityC1: 2.7255E-03
+sutherlandConductivityS: 2.7027E+02
\ No newline at end of file
diff --git a/fluidDatabase/O.dat b/fluidDatabase/O.dat
new file mode 100644
index 0000000..cfacf50
--- /dev/null
+++ b/fluidDatabase/O.dat
@@ -0,0 +1,17 @@
+# general fluid properties
+n: 1.5
+molarMass: 15.9994
+vibrationalTemperature: []
+
+# Sutherland's coefficients for transport properties
+# ----- viscosity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandViscosityC1: 1.9664E-06
+sutherlandViscosityS: 1.1649E+02
+
+# ----- thermal conductivity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandConductivityC1: 3.8319E-03
+sutherlandConductivityS: 1.1649E+02
\ No newline at end of file
diff --git a/fluidDatabase/O2.dat b/fluidDatabase/O2.dat
new file mode 100644
index 0000000..46a789d
--- /dev/null
+++ b/fluidDatabase/O2.dat
@@ -0,0 +1,17 @@
+# general fluid properties
+n: 2.5
+molarMass: 31.9988
+vibrationalTemperature: [2273]
+
+# Sutherland's coefficients for transport properties
+# ----- viscosity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandViscosityC1: 1.7146E-06
+sutherlandViscosityS: 1.3610E+02
+
+# ----- thermal conductivity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandConductivityC1: 3.0048E-03
+sutherlandConductivityS: 3.0610E+02
\ No newline at end of file
diff --git a/fluidDatabase/OH.dat b/fluidDatabase/OH.dat
new file mode 100644
index 0000000..0560df1
--- /dev/null
+++ b/fluidDatabase/OH.dat
@@ -0,0 +1,17 @@
+# general fluid properties
+n: 2.5
+molarMass: 17.0073
+vibrationalTemperature: [5374.2]
+
+# Sutherland's coefficients for transport properties
+# ----- viscosity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandViscosityC1: 2.0274E-06
+sutherlandViscosityS: 1.1649E+02
+
+# ----- thermal conductivity -----
+# values obtained from least squares curve fit of data
+# from [250.0 K - 1000.0 K]
+sutherlandConductivityC1: 4.8939E-03
+sutherlandConductivityS: 1.4471E+02
\ No newline at end of file
diff --git a/fluidDatabase/air.dat b/fluidDatabase/air.dat
new file mode 100644
index 0000000..e3cd205
--- /dev/null
+++ b/fluidDatabase/air.dat
@@ -0,0 +1,12 @@
+# general fluid properties
+n: 2.5
+molarMass: 28.97
+vibrationalTemperature: [3056.0]
+
+# transport properties -- viscosity
+sutherlandViscosityC1: 1.458e-6
+sutherlandViscosityS: 110.4
+
+# transport properties -- thermal conductivity
+sutherlandConductivityC1: 2.495e-3
+sutherlandConductivityS: 194.0
\ No newline at end of file
diff --git a/include/arrayView.hpp b/include/arrayView.hpp
new file mode 100644
index 0000000..b47bf3d
--- /dev/null
+++ b/include/arrayView.hpp
@@ -0,0 +1,451 @@
+/* This file is part of aither.
+ Copyright (C) 2015-18 Michael Nucci (michael.nucci@gmail.com)
+
+ Aither is free software: you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ Aither is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with this program. If not, see . */
+
+#ifndef ARRAYVIEWHEADERDEF
+#define ARRAYVIEWHEADERDEF
+
+#include
+#include
+#include
+#include
+#include
+#include
+#include
+#include
+#include "macros.hpp"
+#include "varArray.hpp"
+#include "vector3d.hpp"
+#include "physicsModels.hpp"
+#include "conserved.hpp"
+
+using std::ostream;
+using std::vector;
+using std::endl;
+using std::unique_ptr;
+
+/* class to store a view of an array. This is useful to slice out data from a
+ * std::vector.
+ */
+
+// forward class declarations
+class primitive;
+class residual;
+
+template
+class arrayView {
+ static_assert(std::is_base_of::value,
+ "arrayView requires T1 to be varArray type!");
+ static_assert(std::is_arithmetic::value,
+ "arrayView requires T2 to be an arithmetic type!");
+
+ typename vector::const_iterator begin_;
+ typename vector::const_iterator end_;
+ int momentumIndex_;
+ int energyIndex_;
+ int turbulenceIndex_;
+
+ public:
+ // constructor
+ arrayView(const typename vector::const_iterator &b,
+ const typename vector::const_iterator &e, const int &numSpecies)
+ : begin_(b),
+ end_(e),
+ momentumIndex_(numSpecies),
+ energyIndex_(momentumIndex_ + 3),
+ turbulenceIndex_(energyIndex_ + 1) {}
+
+ // move constructor and assignment operator
+ arrayView(arrayView&&) noexcept = default;
+ arrayView& operator=(arrayView&&) noexcept = default;
+
+ // copy constructor and assignment operator
+ arrayView(const arrayView &) = default;
+ arrayView &operator=(const arrayView &) = default;
+
+ // member functions
+ T2 Sum() const { return std::accumulate(begin_, end_, T2(0)); }
+ T1 CopyData() const { return {begin_, end_, this->NumSpecies()}; }
+ T1 GetViewType() const { return T1(this->Size(), this->NumSpecies(), 0.0); }
+ T1 Squared() const { return (*this) * (*this); }
+ auto Size() const { return std::distance(begin_, end_); }
+ auto begin() const { return begin_; }
+ auto end() const { return end_; }
+ int NumSpecies() const { return momentumIndex_; }
+ int NumTurbulence() const { return this->Size() - turbulenceIndex_; }
+ bool IsMultiSpecies() const { return this->NumSpecies() > 1; }
+ bool HasTurbulenceData() const { return this->Size() != turbulenceIndex_; }
+ int MomentumXIndex() const { return momentumIndex_; }
+ int MomentumYIndex() const { return momentumIndex_ + 1; }
+ int MomentumZIndex() const { return momentumIndex_ + 2; }
+ int EnergyIndex() const { return energyIndex_; }
+ int TurbulenceIndex() const { return turbulenceIndex_; }
+ T2 SpeciesSum() const {
+ return std::accumulate(begin_, begin_ + this->NumSpecies(), T2(0));
+ }
+ const T2 &SpeciesN(const int &ii) const {
+ MSG_ASSERT(ii < momentumIndex_, "requesting species variable out of range");
+ return (*this)[ii];
+ }
+ const T2 &MomentumX() const { return (*this)[momentumIndex_]; }
+ const T2 &MomentumY() const { return (*this)[momentumIndex_ + 1]; }
+ const T2 &MomentumZ() const { return (*this)[momentumIndex_ + 2]; }
+ const T2 &Energy() const { return (*this)[energyIndex_]; }
+ const T2 &TurbulenceN(const int &ii) const {
+ MSG_ASSERT(turbulenceIndex_ + ii < this->Size(),
+ "requesting turbulence variable out of range");
+ return (*this)[turbulenceIndex_ + ii];
+ }
+
+ arrayView GetView() const { return *this; }
+
+ // --------------------------------------------------------------------------
+ // getters for primitives ---------------------------------------------------
+ const T2 &RhoN(const int &ii) const {
+ static_assert(std::is_same::value ||
+ std::is_same::value,
+ "getter only valid for primitive/conserved type!");
+ return this->SpeciesN(ii);
+ }
+ T2 Rho() const {
+ static_assert(std::is_same::value ||
+ std::is_same::value,
+ "getter only valid for primitive/conserved type!");
+ return this->SpeciesSum();
+ }
+ T2 MassFractionN(const int &ii) const {
+ static_assert(std::is_same::value ||
+ std::is_same::value,
+ "getter only valid for primitive/conserved type!");
+ return this->RhoN(ii) / this->Rho();
+ }
+ vector RhoVec() const {
+ static_assert(std::is_same::value ||
+ std::is_same::value,
+ "getter only valid for primitive/conserved type!");
+ return {this->begin(), this->begin() + this->NumSpecies()};
+ }
+ vector MassFractions() const {
+ static_assert(std::is_same::value ||
+ std::is_same::value,
+ "getter only valid for primitive/conserved type!");
+ vector mf(this->NumSpecies());
+ const auto totalRho = this->Rho();
+ for (auto ii = 0U; ii < mf.size(); ++ii) {
+ mf[ii] = this->RhoN(ii) / totalRho;
+ }
+ return mf;
+ }
+ vector VolumeFractions(const unique_ptr &trans) const {
+ static_assert(std::is_same::value,
+ "getter only valid for primitive type!");
+ // volume fractions equal to mole fractions
+ return trans->MoleFractions(this->MassFractions());
+ }
+ const T2 &U() const {
+ static_assert(std::is_same::value,
+ "getter only valid for primitive type!");
+ return this->MomentumX();
+ }
+ const T2 &V() const {
+ static_assert(std::is_same::value,
+ "getter only valid for primitive type!");
+ return this->MomentumY();
+ }
+ const T2 &W() const {
+ static_assert(std::is_same::value,
+ "getter only valid for primitive type!");
+ return this->MomentumZ();
+ }
+ const T2 &P() const {
+ static_assert(std::is_same::value,
+ "getter only valid for primitive type!");
+ return this->Energy();
+ }
+ const T2 &Tke() const {
+ static_assert(std::is_same::value,
+ "getter only valid for primitive type!");
+ return this->TurbulenceN(0);
+ }
+ const T2 &Omega() const {
+ static_assert(std::is_same::value,
+ "getter only valid for primitive type!");
+ return this->TurbulenceN(1);
+ }
+ const T2 &TurbN(const int &ii) const {
+ static_assert(std::is_same::value,
+ "getter only valid for primitive type!");
+ return this->TurbulenceN(ii);
+ }
+ vector3d Velocity() const {
+ static_assert(std::is_same::value,
+ "getter only valid for primitive type!");
+ return {this->U(), this->V(), this->W()};
+ }
+ T2 SoS(const physics &phys) const;
+ T2 Energy(const physics &phys) const;
+ T2 SpeciesEnthalpy(const physics &phys, const int &ss) const;
+ T2 Enthalpy(const physics &phys) const;
+ T2 Temperature(const unique_ptr &eqnState) const {
+ static_assert(std::is_same::value,
+ "getter only valid for primitive type!");
+ return eqnState->Temperature(this->P(), this->RhoVec());
+ }
+ conserved ConsVars(const physics &phys) const {
+ static_assert(std::is_same::value,
+ "function only valid for primitive type!");
+ return PrimToCons((*this), phys);
+ }
+ T1 UpdateWithConsVars(const physics &phys,
+ const arrayView &du) const {
+ static_assert(std::is_same::value,
+ "function only valid for primitive type!");
+ return UpdatePrimWithCons((*this), phys, du);
+ }
+
+ // --------------------------------------------------------------------------
+ // operator overloads
+ const T2 & operator[](const int &r) const { return *(begin_ + r); }
+
+ inline T1 operator+(const T2 &s) const {
+ auto lhs = this->CopyData();
+ return lhs += s;
+ }
+ inline T1 operator-(const T2 &s) const {
+ auto lhs = this->CopyData();
+ return lhs -= s;
+ }
+ inline T1 operator*(const T2 &s) const {
+ auto lhs = this->CopyData();
+ return lhs *= s;
+ }
+ inline T1 operator/(const T2 &s) const {
+ auto lhs = this->CopyData();
+ return lhs /= s;
+ }
+
+ // destructor
+ ~arrayView() noexcept {}
+};
+
+// function declarations --------------------------------------
+template
+inline const T1 operator+(const arrayView &lhs,
+ const arrayView &rhs) {
+ auto ll = lhs.CopyData();
+ return ll + rhs;
+}
+
+template
+inline const T1 operator-(const arrayView &lhs,
+ const arrayView &rhs) {
+ auto ll = lhs.CopyData();
+ return ll - rhs;
+}
+
+template
+inline const T1 operator*(const arrayView &lhs,
+ const arrayView &rhs) {
+ auto ll = lhs.CopyData();
+ return ll * rhs;
+}
+
+template
+inline const T1 operator/(const arrayView &lhs,
+ const arrayView &rhs) {
+ auto ll = lhs.CopyData();
+ return ll / rhs;
+}
+
+// operator overloads for type T -------------------------------------
+template
+inline const T1 operator+(const T2 &lhs, const arrayView &rhs) {
+ auto result = rhs.CopyData();
+ return result += lhs;
+}
+
+template
+inline const T1 operator-(const T2 &lhs, const arrayView &rhs) {
+ auto result = rhs.CopyData();
+ return lhs - result;
+}
+
+template
+inline const T1 operator*(const T2 &lhs, const arrayView &rhs) {
+ auto result = rhs.CopyData();
+ return result *= lhs;
+}
+
+template
+inline const T1 operator/(const T2 &lhs, const arrayView &rhs) {
+ auto result = rhs.CopyData();
+ return lhs / rhs;
+}
+
+// operator overloads for varArray type T -------------------------------------
+template
+inline const typename std::enable_if_t<
+ std::is_base_of::value, T3>
+operator+(const T3 &lhs, const arrayView &rhs) {
+ T3 result(rhs.begin(), rhs.end(), rhs.NumSpecies());
+ return lhs + result;
+}
+
+template
+inline const typename std::enable_if_t<
+ std::is_base_of::value, T3>
+operator-(const T3 &lhs, const arrayView &rhs) {
+ T3 result(rhs.begin(), rhs.end(), rhs.NumSpecies());
+ return lhs - result;
+}
+
+template
+inline const typename std::enable_if_t<
+ std::is_base_of::value, T3>
+operator*(const T3 &lhs, const arrayView &rhs) {
+ T3 result(rhs.begin(), rhs.end(), rhs.NumSpecies());
+ return lhs * result;
+}
+
+template
+inline const typename std::enable_if_t<
+ std::is_base_of::value, T3>
+operator/(const T3 &lhs, const arrayView &rhs) {
+ T3 result(rhs.begin(), rhs.end(), rhs.NumSpecies());
+ return lhs / result;
+}
+
+// ---------------------------------------------------------------------------
+template
+inline const typename std::enable_if_t<
+ std::is_base_of::value, T3>
+operator+(const arrayView &lhs, const T3 &rhs) {
+ T3 result(lhs.begin(), lhs.end(), lhs.NumSpecies());
+ return result += rhs;
+}
+
+template
+inline const typename std::enable_if_t<
+ std::is_base_of::value, T3>
+operator-(const arrayView &lhs, const T3 &rhs) {
+ T3 result(lhs.begin(), lhs.end(), lhs.NumSpecies());
+ return result -= rhs;
+}
+
+template
+inline const typename std::enable_if_t<
+ std::is_base_of::value, T3>
+operator*(const arrayView &lhs, const T3 &rhs) {
+ T3 result(lhs.begin(), lhs.end(), lhs.NumSpecies());
+ return result *= rhs;
+}
+
+template
+inline const typename std::enable_if_t<
+ std::is_base_of::value, T3>
+operator/(const arrayView &lhs, const T3 &rhs) {
+ T3 result(lhs.begin(), lhs.end(), lhs.NumSpecies());
+ return result /= rhs;
+}
+
+
+// ---------------------------------------------------------------------------
+// operation overload for << - allows use of cout, cerr, etc.
+template
+ostream &operator<<(ostream &os, const arrayView &m) {
+ for (auto rr = 0; rr < m.Size(); rr++) {
+ os << m[rr] << endl;
+ }
+ return os;
+}
+
+// using typedefs
+using varArrayView = arrayView;
+using primitiveView = arrayView;
+using conservedView = arrayView;
+using residualView = arrayView;
+
+
+template
+double SpeedOfSound(const T &state, const physics &phys) {
+ static_assert(std::is_same::value ||
+ std::is_same::value,
+ "T requires primitive or primativeView type");
+ return sqrt(phys.Thermodynamic()->Gamma(state.Temperature(phys.EoS()),
+ state.MassFractions()) *
+ state.P() / state.Rho());
+}
+
+template
+T2 arrayView::SoS(const physics &phys) const {
+ static_assert(std::is_same::value,
+ "getter only valid for primitive type!");
+ return SpeedOfSound(*this, phys);
+}
+
+template
+double EnthalpyFunc(const T &state, const physics &phys) {
+ static_assert(std::is_same::value ||
+ std::is_same::value,
+ "T requires primitive or primativeView type");
+ const auto t = state.Temperature(phys.EoS());
+ return phys.EoS()->Enthalpy(phys.Thermodynamic(), t, state.Velocity().Mag(),
+ state.MassFractions());
+}
+
+template
+double SpeciesEnthalpyFunc(const T &state, const physics &phys, const int &ss) {
+ static_assert(std::is_same::value ||
+ std::is_same::value,
+ "T requires primitive or primativeView type");
+ const auto t = state.Temperature(phys.EoS());
+ return phys.EoS()->SpeciesEnthalpy(phys.Thermodynamic(), t,
+ state.Velocity().Mag(), ss);
+}
+
+template
+T2 arrayView::Enthalpy(const physics &phys) const {
+ static_assert(std::is_same::value,
+ "getter only valid for primitive type!");
+ return EnthalpyFunc(*this, phys);
+}
+
+template
+T2 arrayView::SpeciesEnthalpy(const physics &phys,
+ const int &ss) const {
+ static_assert(std::is_same::value,
+ "getter only valid for primitive type!");
+ return SpeciesEnthalpyFunc(*this, phys, ss);
+}
+
+template
+double InternalEnergy(const T &state, const physics &phys) {
+ static_assert(std::is_same::value ||
+ std::is_same::value,
+ "T requires primitive or primativeView type");
+ const auto t = state.Temperature(phys.EoS());
+ return phys.EoS()->Energy(
+ phys.EoS()->SpecEnergy(phys.Thermodynamic(), t, state.MassFractions()),
+ state.Velocity().Mag());
+}
+
+template
+T2 arrayView::Energy(const physics &phys) const {
+ static_assert(std::is_same::value,
+ "getter only valid for primitive type!");
+ return InternalEnergy(*this, phys);
+}
+
+#endif
diff --git a/include/blkMultiArray3d.hpp b/include/blkMultiArray3d.hpp
new file mode 100644
index 0000000..012d4e8
--- /dev/null
+++ b/include/blkMultiArray3d.hpp
@@ -0,0 +1,432 @@
+/* This file is part of aither.
+ Copyright (C) 2015-18 Michael Nucci (michael.nucci@gmail.com)
+
+ Aither is free software: you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ Aither is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with this program. If not, see . */
+
+#ifndef BLKMULTIARRAY3DHEADERDEF
+#define BLKMULTIARRAY3DHEADERDEF
+
+/* This file contains the header and implementation for a multidimensional (3D)
+ array class. The class is to act as a container to store all data types,
+ and provide easy access to elements with i, j, k indexing.
+ */
+
+#include // ostream
+#include // vector
+#include // string
+#include // unique_ptr
+#include "mpi.h"
+#include "vector3d.hpp"
+#include "boundaryConditions.hpp" // connection
+#include "range.hpp" // range
+#include "multiArray3d.hpp"
+#include "varArray.hpp"
+#include "arrayView.hpp"
+
+using std::ostream;
+using std::endl;
+using std::cout;
+using std::cerr;
+using std::vector;
+using std::string;
+using std::unique_ptr;
+
+template
+class blkMultiArray3d : public multiArray3d {
+ static_assert(std::is_base_of::value,
+ "blkMultiArray3d requires a varArray based type");
+
+ int numSpecies_;
+
+ public:
+ // constructor
+ blkMultiArray3d(const int &ii, const int &jj, const int &kk, const int &ng,
+ const T &init)
+ : multiArray3d(ii, jj, kk, ng, init.Size()),
+ numSpecies_(init.NumSpecies()) {
+ for (auto kk = this->StartK(); kk < this->EndK(); kk++) {
+ for (auto jj = this->StartJ(); jj < this->EndJ(); jj++) {
+ for (auto ii = this->StartI(); ii < this->EndI(); ii++) {
+ for (auto bb = 0; bb < this->BlockSize(); ++bb) {
+ (*this)(ii, jj, kk, bb) = init[bb];
+ }
+ }
+ }
+ }
+ }
+ blkMultiArray3d(const int &ii, const int &jj, const int &kk, const int &ng,
+ const int &bs, const int &ns, const double &init)
+ : multiArray3d(ii, jj, kk, ng, bs, init), numSpecies_(ns) {}
+ blkMultiArray3d(const int &ii, const int &jj, const int &kk, const int &ng,
+ const int &bs, const int &ns)
+ : multiArray3d(ii, jj, kk, ng, bs), numSpecies_(ns) {}
+ blkMultiArray3d(const int &ii, const int &jj, const int &kk, const int &ng,
+ const std::pair &info)
+ : multiArray3d(ii, jj, kk, ng, info),
+ numSpecies_(info.second) {}
+ blkMultiArray3d() : multiArray3d(), numSpecies_(0) {}
+
+ // move constructor and assignment operator
+ blkMultiArray3d(blkMultiArray3d &&) noexcept = default;
+ blkMultiArray3d &operator=(blkMultiArray3d &&) noexcept = default;
+
+ // copy constructor and assignment operator
+ blkMultiArray3d(const blkMultiArray3d &) = default;
+ blkMultiArray3d &operator=(const blkMultiArray3d &) = default;
+
+ // member functions
+ std::pair BlockInfo() const override {
+ return std::make_pair(this->BlockSize(), numSpecies_);
+ }
+
+ template
+ void InsertBlock(const int &ii, const int &jj, const int &kk, const T2 &arr) {
+ static_assert(std::is_same::value ||
+ std::is_same, T2>::value,
+ "blkMultiArray3d requires a varArray based type!");
+ MSG_ASSERT(arr.Size() == this->BlockSize(), "block size must match");
+ auto start = this->GetLoc1D(ii, jj, kk);
+ std::copy(arr.begin(), arr.end(), this->begin() + start);
+ }
+
+ template
+ void InsertBlock(const string &dir, const int &d1, const int &d2,
+ const int &d3, const T2 &arr) {
+ static_assert(std::is_same::value ||
+ std::is_same, T2>::value,
+ "blkMultiArray3d requires a varArray based type!");
+ MSG_ASSERT(arr.Size() == this->BlockSize(), "block size must match");
+
+ auto start = 0;
+ if (dir == "i") { // direction 1 is i
+ start = this->GetLoc1D(d1, d2, d3);
+ } else if (dir == "j") { // direction 1 is j
+ start = this->GetLoc1D(d3, d1, d2);
+ } else if (dir == "k") { // direction 1 is k
+ start = this->GetLoc1D(d2, d3, d1);
+ } else {
+ cerr << "ERROR: Direction " << dir << " is not recognized!" << endl;
+ exit(EXIT_FAILURE);
+ }
+ std::copy(arr.begin(), arr.end(), this->begin() + start);
+ }
+
+ blkMultiArray3d Slice(const range &, const range &, const range &) const;
+ blkMultiArray3d Slice(const string &, const range &,
+ const bool = false) const;
+ blkMultiArray3d Slice(const string &, int, int, const bool = false,
+ const string = "cell", const bool = false,
+ const bool = false) const;
+ blkMultiArray3d Slice(const string &, int, range, range,
+ const string = "cell", const int = 0) const;
+
+ template
+ void Insert(const range &, const range &, const range &, const TT &);
+ template
+ void Insert(const string &, const range &, const TT &, const bool = false);
+ template
+ void Insert(const string &, int, int, const TT &, const bool = false,
+ const string = "cell", const bool = false, const bool = false);
+ template
+ void Insert(const string &, int, range, range, const TT &,
+ const string = "cell", const int = 0);
+
+ template
+ void SwapSlice(const connection &conn, TT &array);
+ void SwapSliceMPI(const connection &conn, const int &rank,
+ const MPI_Datatype &MPI_arrData, const int tag = 1);
+ template
+ void PutSlice(const TT &, const connection &, const int &);
+
+ T GetCopy(const int &ii, const int &jj, const int &kk) const {
+ auto start = this->GetLoc1D(ii, jj, kk);
+ return {this->begin() + start, this->begin() + start + this->BlockSize(),
+ numSpecies_};
+ }
+
+ // operator overloads
+ arrayView operator()(const int &ii, const int &jj,
+ const int &kk) const {
+ auto start = this->GetLoc1D(ii, jj, kk);
+ return {this->begin() + start, this->begin() + start + this->BlockSize(),
+ numSpecies_};
+ }
+ arrayView operator()(const int &ind) const {
+ auto start = ind * this->BlockSize();
+ return {this->begin() + start, this->begin() + start + this->BlockSize(),
+ numSpecies_};
+ }
+ arrayView operator()(const string &dir, const int &d1,
+ const int &d2, const int &d3) const {
+ auto start = 0;
+ if (dir == "i") { // direction 1 is i
+ start = this->GetLoc1D(d1, d2, d3);
+ } else if (dir == "j") { // direction 1 is j
+ start = this->GetLoc1D(d3, d1, d2);
+ } else if (dir == "k") { // direction 1 is k
+ start = this->GetLoc1D(d2, d3, d1);
+ } else {
+ cerr << "ERROR: Direction " << dir << " is not recognized!" << endl;
+ exit(EXIT_FAILURE);
+ }
+ return {this->begin() + start, this->begin() + start + this->BlockSize(),
+ numSpecies_};
+ }
+ double &operator()(const int &ii, const int &jj, const int &kk,
+ const int &bb) {
+ MSG_ASSERT(bb <= this->BlockSize(), "accessing data out of block index");
+ return (*this)[this->GetLoc1D(ii, jj, kk) + bb];
+ }
+ const double &operator()(const int &ii, const int &jj, const int &kk,
+ const int &bb) const {
+ MSG_ASSERT(bb <= this->BlockSize(), "accessing data out of block index");
+ return (*this)[this->GetLoc1D(ii, jj, kk) + bb];
+ }
+
+ void ClearResize(const int &ii, const int &jj, const int &kk,
+ const int &ng, const int &bs, const int &ns) {
+ *this = blkMultiArray3d(ii, jj, kk, ng, bs, ns);
+ }
+ void ClearResize(const int &ii, const int &jj, const int &kk, const int &ng,
+ const int &bs, const int &ns, const T &val) {
+ *this = blkMultiArray3d(ii, jj, kk, ng, bs, ns, val);
+ }
+
+ blkMultiArray3d GrowI() const;
+ blkMultiArray3d GrowJ() const;
+ blkMultiArray3d GrowK() const;
+
+ // destructor
+ ~blkMultiArray3d() noexcept {}
+};
+
+// ---------------------------------------------------------------------------
+// member function definitions
+
+// operation overload for << - allows use of cout, cerr, etc.
+template
+ostream &operator<<(ostream &os, const blkMultiArray3d &arr) {
+ os << "Size: " << arr.NumI() << ", " << arr.NumJ() << ", " << arr.NumK()
+ << endl;
+ os << "Number of ghost layers: " << arr.GhostLayers() << endl;
+
+ for (auto kk = arr.StartK(); kk < arr.EndK(); ++kk) {
+ for (auto jj = arr.StartJ(); jj < arr.EndJ(); ++jj) {
+ for (auto ii = arr.StartI(); ii < arr.EndI(); ++ii) {
+ os << ii << ", " << jj << ", " << kk << ": " << arr(ii, jj, kk) << endl;
+ }
+ }
+ }
+ return os;
+}
+
+// member function to return a slice of the array
+// this is the main slice function that all other overloaded slice functions call
+template
+blkMultiArray3d blkMultiArray3d::Slice(const range &ir, const range &jr,
+ const range &kr) const {
+ // ir -- i-index range to take slice [inclusive, exclusive)
+ // jr -- j-index range to take slice [inclusive, exclusive)
+ // kr -- k-index range to take slice [inclusive, exclusive)
+ return SliceArray((*this), ir, jr, kr);
+}
+
+// member function to return a slice of the array
+// Overload to slice only in one direction. Given a 3D array, this slice returns
+// a plane with normal direction dir, or a smaller 3D array where the direction
+// dir is sliced over dirRange. It also has the ability to include or ignore
+// ghost cells in its planar slices
+template
+blkMultiArray3d blkMultiArray3d::Slice(const string &dir,
+ const range &dirRange,
+ const bool physOnly) const {
+ // dir -- direction of slice
+ // dirRange -- range of slice in direction given
+ // phsOnly -- flag to only include physical cells in the two directions that
+ // are not specified as dir
+ return SliceArray((*this), dir, dirRange, physOnly);
+}
+
+// member function to return a slice of the array
+// overload to slice line out of array
+template
+blkMultiArray3d blkMultiArray3d::Slice(const string &dir, int d2Ind,
+ int d3Ind, const bool physOnly,
+ const string id, const bool upper2,
+ const bool upper3) const {
+ // dir -- direction of line slice (direction 1)
+ // d2Ind -- index of direction 2
+ // d3Ind -- index of direction 3
+ // physOnly -- flag to only include physical cells in line slice
+ // id -- type of multiArray3d being sliced: cell, i, j, or k
+ // d2Ind and d3Ind are supplied as cell indices, but may need to be
+ // altered if the array is storing i, j, or k face data
+ // upper2 -- flag to determine if direction 2 is at upper index
+ // upper3 -- flag to determine if direction 3 is at upper index
+ return SliceArray((*this), dir, d2Ind, d3Ind, physOnly, id, upper2, upper3);
+}
+
+// overload to slice plane out of array
+// Identical to previous slice overload, but more general in that in can slice
+// over a subset of direction 2 & 3. This is useful to slice out a plane that
+// borders a boundary condition patch.
+template
+blkMultiArray3d blkMultiArray3d::Slice(const string &dir, int dirInd,
+ range dir1, range dir2,
+ const string id,
+ const int type) const {
+ // dir -- normal direction of planar slice
+ // dirInd -- index in normal direction
+ // dir1 -- range of direction 1 (direction 3 is normal to slice)
+ // dir2 -- range of direction 2 (direction 3 is normal to slice)
+ // id -- id of array being sliced (i, j, k for faces, cell for cells)
+ // type -- surface type of dir
+ return SliceArray((*this), dir, dirInd, dir1, dir2, id, type);
+}
+
+// member function to insert an array into this one
+// this is the main insert funciton that all other overloaded insert functions
+// call
+template
+template
+void blkMultiArray3d::Insert(const range &ir, const range &jr,
+ const range &kr, const TT &arr) {
+ // ir -- i-index range to take slice [inclusive, exclusive)
+ // jr -- j-index range to take slice [inclusive, exclusive)
+ // kr -- k-index range to take slice [inclusive, exclusive)
+ // arr -- array to insert into this one
+ InsertArray((*this), ir, jr, kr, arr);
+}
+
+// Overload to insert only in one direction. Given a 3D array, this inserts a
+// plane with normal direction dir, or a smaller 3D array where the direction
+// dir is inserted over dirRange. It also has the ability to include or ignore
+// ghost cells in its planar inserts
+template
+template
+void blkMultiArray3d::Insert(const string &dir, const range &dirRange,
+ const TT &arr, const bool physOnly) {
+ // dir -- direction of slice to insert
+ // dirRange -- range to insert slice into in direction given
+ // arr -- array to insert
+ // phsOnly -- flag to only include physical cells in the two directions that
+ // are not specified as dir
+ InsertArray((*this), dir, dirRange, arr, physOnly);
+}
+
+// overload to insert line into array
+template
+template
+void blkMultiArray3d::Insert(const string &dir, int d2Ind, int d3Ind,
+ const TT &arr, const bool physOnly,
+ const string id, const bool upper2,
+ const bool upper3) {
+ // dir -- direction of line slice to insert (direction 1)
+ // d2Ind -- index of direction 2 to insert into
+ // d3Ind -- index of direction 3 to insert into
+ // physOnly -- flag to only include physical cells in line insert
+ // id -- type of multiArray3d being sliced: cell, i, j, or k
+ // d2Ind and d3Ind are supplied as cell indices, but may need to be
+ // altered if the array is storing i, j, or k face data
+ // upper2 -- flag to determine if direction 2 is at upper index
+ // upper3 -- flag to determine if direction 3 is at upper index
+ InsertArray((*this), dir, d2Ind, d3Ind, arr, physOnly, id, upper2, upper3);
+}
+
+// overload to insert plane into array
+// Identical to previous insert overload, but more general in that in can insert
+// over a subset of direction 2 & 3. This is useful to insert into a plane that
+// borders a boundary condition patch.
+template
+template
+void blkMultiArray3d::Insert(const string &dir, int dirInd, range dir1,
+ range dir2, const TT &arr, const string id,
+ const int type) {
+ // dir -- normal direction of planar slice
+ // dirInd -- index in normal direction
+ // dir1 -- range of direction 1 (direction 3 is normal to slice)
+ // dir2 -- range of direction 2 (direction 3 is normal to slice)
+ // arr -- array to insert
+ // id -- id of array being inserted into (i, j, k for faces, cell for cells)
+ // type -- surface type of dir
+ InsertArray((*this), dir, dirInd, dir1, dir2, arr, id, type);
+}
+
+template
+blkMultiArray3d blkMultiArray3d::GrowI() const {
+ return GrowInI(*this);
+}
+template
+blkMultiArray3d blkMultiArray3d::GrowJ() const {
+ return GrowInJ(*this);
+}
+template
+blkMultiArray3d blkMultiArray3d::GrowK() const {
+ return GrowInK(*this);
+}
+
+/* Function to swap slice using MPI. This is similar to the SwapSlice
+ function, but is called when the neighboring procBlocks are on different
+ processors.
+*/
+template
+void blkMultiArray3d::SwapSliceMPI(const connection &conn, const int &rank,
+ const MPI_Datatype &MPI_arrData,
+ const int tag) {
+ // conn -- connection boundary information
+ // rank -- processor rank
+ // MPI_arrData -- MPI datatype for passing data in *this
+ // tag -- id for MPI swap (default 1)
+ SwapSliceParallel((*this), conn, rank, MPI_arrData, tag);
+}
+
+/* Function to swap ghost cells between two blocks at an connection
+boundary. Slices are removed from the physical cells (extending into ghost cells
+at the edges) of one block and inserted into the ghost cells of its partner
+block. The reverse is also true. The slices are taken in the coordinate system
+orientation of their parent block.
+
+ Interior Cells Ghost Cells Ghost Cells Interior Cells
+ ________ ______|________ _________ _______________|_______ _________
+Ui-3/2 Ui-1/2 | Uj+1/2 Uj+3/2 Ui-3/2 Ui-1/2 | Uj+1/2 Uj+3/2
+ | | | | | | | | | |
+ | Ui-1 | Ui | Uj | Uj+1 | | Ui-1 | Ui | Uj | Uj+1 |
+ | | | | | | | | | |
+ |________|______|________|_________| |________|______|_______|_________|
+ | |
+
+The above diagram shows the resulting values after the ghost cell swap. The
+logic ensures that the ghost cells at the connection boundary exactly match
+their partner block as if there were no separation in the grid.
+*/
+template
+template
+void blkMultiArray3d::SwapSlice(const connection &conn, TT &array) {
+ // conn -- connection boundary information
+ // array -- second array involved in connection boundary
+ SwapSliceLocal((*this), conn, array);
+}
+
+template
+template
+void blkMultiArray3d::PutSlice(const TT &array, const connection &inter,
+ const int &d3) {
+ // array -- array to insert into *this
+ // inter -- connection data structure defining patches and orientation
+ // d3 -- distance of direction normal to patch to insert
+ InsertSlice((*this), array, inter, d3);
+}
+
+
+#endif
diff --git a/include/boundaryConditions.hpp b/include/boundaryConditions.hpp
index b946e98..4b0fb0f 100644
--- a/include/boundaryConditions.hpp
+++ b/include/boundaryConditions.hpp
@@ -1,5 +1,5 @@
/* This file is part of aither.
- Copyright (C) 2015-17 Michael Nucci (michael.nucci@gmail.com)
+ Copyright (C) 2015-18 Michael Nucci (michael.nucci@gmail.com)
Aither is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
@@ -29,6 +29,9 @@ one block. */
#include // string
#include // ostream
#include // unique_ptr
+#include