diff --git a/examples/inviscid-flow-sw/UpwindContinuity.h b/examples/inviscid-flow-sw/UpwindContinuity.h
index a375f8b9597f735873a9a148e6fb3a3aa9244830..19150164f063d8df12f94a4ea9b26ed591bf3487 100644
--- a/examples/inviscid-flow-sw/UpwindContinuity.h
+++ b/examples/inviscid-flow-sw/UpwindContinuity.h
@@ -171,10 +171,6 @@ class UpwindContinuity< Meshes::Grid< 1, MeshReal, Device, MeshIndex >, Real, In
- this->positiveDensityFlux( u[ west ], velocity_x_west , pressure_west )
- this->negativeDensityFlux( u[ center ], velocity_x_center, pressure_center )
+ this->negativeDensityFlux( u[ east ], velocity_x_east , pressure_east )
-// u[ center ] / ( 2 * this->gamma ) * ( ( 2 * this->gamma - 1 ) * velocity_x_center + std::sqrt( this->multiply( this->gamma, pressure_center ) / u[ center ] ) )
-// - u[ west ] / ( 2 * this->gamma ) * ( ( 2 * this->gamma - 1 ) * velocity_x_west + std::sqrt( this->gamma * pressure_west / u[ west ] ) )
-// - u[ center ] / ( 2 * this->gamma ) * ( velocity_x_center - std::sqrt( this->gamma * pressure_center / u[ center ] ) )
-// + u[ east ] / ( 2 * this->gamma ) * ( velocity_x_east - std::sqrt( this->gamma * pressure_east / u[ east ] ) )
);
}
@@ -251,17 +247,17 @@ class UpwindContinuity< Meshes::Grid< 2, MeshReal, Device, MeshIndex >, Real, In
const RealType& velocity_y_north = this->velocity.template getData< DeviceType >()[ 1 ].template getData< DeviceType >()[ north ];
const RealType& velocity_y_south = this->velocity.template getData< DeviceType >()[ 1 ].template getData< DeviceType >()[ south ];
- return -hxInverse * (
- u[ center ] / ( 2 * this->gamma ) * ( ( 2 * this->gamma - 1 ) * velocity_x_center + std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- - u[ west ] / ( 2 * this->gamma ) * ( ( 2 * this->gamma - 1 ) * velocity_x_west + std::sqrt( this->gamma * pressure_west / u[ west ] ) )
- - u[ center ] / ( 2 * this->gamma ) * ( velocity_x_center - std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- + u[ east ] / ( 2 * this->gamma ) * ( velocity_x_east - std::sqrt( this->gamma * pressure_east / u[ east ] ) )
+ return -hxInverse * (
+ this->positiveDensityFlux( u[ center ], velocity_x_center, pressure_center )
+ - this->positiveDensityFlux( u[ west ], velocity_x_west , pressure_west )
+ - this->negativeDensityFlux( u[ center ], velocity_x_center, pressure_center )
+ + this->negativeDensityFlux( u[ east ], velocity_x_east , pressure_east )
)
- -hyInverse * (
- u[ center ] / ( 2 * this->gamma ) * ( ( 2 * this->gamma - 1 ) * velocity_y_center + std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- - u[ south ] / ( 2 * this->gamma ) * ( ( 2 * this->gamma - 1 ) * velocity_y_south + std::sqrt( this->gamma * pressure_south / u[ south ] ) )
- - u[ center ] / ( 2 * this->gamma ) * ( velocity_y_center - std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- + u[ north ] / ( 2 * this->gamma ) * ( velocity_y_north - std::sqrt( this->gamma * pressure_north / u[ north ] ) )
+ -hyInverse * (
+ this->positiveDensityFlux( u[ center ], velocity_y_center, pressure_center )
+ - this->positiveDensityFlux( u[ south ], velocity_y_south , pressure_south )
+ - this->negativeDensityFlux( u[ center ], velocity_y_center, pressure_center )
+ + this->negativeDensityFlux( u[ north ], velocity_y_north , pressure_north )
);
}
diff --git a/examples/inviscid-flow-sw/UpwindEnergy.h b/examples/inviscid-flow-sw/UpwindEnergy.h
index 29d74e32e94c46e2e5760195a4c910731985cf18..44792f8b3c7ae8094b4549a47c939098877dd5d7 100644
--- a/examples/inviscid-flow-sw/UpwindEnergy.h
+++ b/examples/inviscid-flow-sw/UpwindEnergy.h
@@ -271,89 +271,17 @@ class UpwindEnergy< Meshes::Grid< 2, MeshReal, Device, MeshIndex >, Real, Index
const RealType& velocity_y_north = this->velocity.template getData< DeviceType >()[ 1 ].template getData< DeviceType >()[ north ];
const RealType& velocity_y_south = this->velocity.template getData< DeviceType >()[ 1 ].template getData< DeviceType >()[ south ];
- return -hxInverse * (
- ( u[ center ] / ( 2 * this->gamma ) )
- * (
- ( ( 2 * this->gamma - 1 ) * velocity_x_center + std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- * ( velocity_x_center * velocity_x_center + velocity_y_center * velocity_y_center ) / 2
- + velocity_x_center * std::sqrt( this->gamma * pressure_center / u[ center ] )
- * ( velocity_x_center + std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- + this->gamma * pressure_center / u[ center ]
- * ( velocity_x_center + std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- / ( this->gamma - 1 )
- )
- - ( u[ west ] / ( 2 * this->gamma ) )
- * (
- ( ( 2 * this->gamma - 1 ) * velocity_x_west + std::sqrt( this->gamma * pressure_west / u[ west ] ) )
- * ( velocity_x_west * velocity_x_west + velocity_y_west * velocity_y_west ) / 2
- + velocity_x_west * std::sqrt( this->gamma * pressure_west / u[ west ] )
- * ( velocity_x_west + std::sqrt( this->gamma * pressure_west / u[ west ] ) )
- + this->gamma * pressure_west / u[ west ]
- * ( velocity_x_west + std::sqrt( this->gamma * pressure_west / u[ west ] ) )
- / ( this->gamma - 1 )
- )
- - ( u[ center ] / ( 2 * this->gamma ) )
- * (
- ( velocity_x_center - std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- * ( velocity_x_center * velocity_x_center + velocity_y_center * velocity_y_center ) / 2
- + velocity_x_center * std::sqrt( this->gamma * pressure_center / u[ center ] )
- * ( std::sqrt( this->gamma * pressure_center / u[ center ] ) - velocity_x_center )
- + this->gamma * pressure_center / u[ center ]
- * ( velocity_x_center - std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- / ( this->gamma - 1 )
- )
- + ( u[ east ] / ( 2 * this->gamma ) )
- * (
- ( velocity_x_east - std::sqrt( this->gamma * pressure_east / u[ east ] ) )
- * ( velocity_x_east * velocity_x_east + velocity_y_east * velocity_y_east ) / 2
- + velocity_x_east * std::sqrt( this->gamma * pressure_east / u[ east ] )
- * ( std::sqrt( this->gamma * pressure_east / u[ east ] ) - velocity_x_east )
- + this->gamma * pressure_east / u[ east ]
- * ( velocity_x_east - std::sqrt( this->gamma * pressure_east / u[ east ] ) )
- / ( this->gamma - 1 )
- )
+ return -hxInverse * (
+ this->positiveEnergyFlux( u[ center ], velocity_x_center, velocity_y_center, 0, pressure_center)
+ - this->positiveEnergyFlux( u[ west ], velocity_x_west , velocity_y_west , 0, pressure_west )
+ - this->negativeEnergyFlux( u[ center ], velocity_x_center, velocity_y_center, 0, pressure_center)
+ + this->negativeEnergyFlux( u[ east ], velocity_x_east , velocity_y_east , 0, pressure_east )
)
- -hyInverse * (
- ( u[ center ] / ( 2 * this->gamma ) )
- * (
- ( ( 2 * this->gamma - 1 ) * velocity_y_center + std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- * ( velocity_x_center * velocity_x_center + velocity_y_center * velocity_y_center ) / 2
- + velocity_y_center * std::sqrt( this->gamma * pressure_center / u[ center ] )
- * ( velocity_y_center + std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- + this->gamma * pressure_center / u[ center ]
- * ( velocity_y_center + std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- / ( this->gamma - 1 )
- )
- - ( u[ south ] / ( 2 * this->gamma ) )
- * (
- ( ( 2 * this->gamma - 1 ) * velocity_y_south + std::sqrt( this->gamma * pressure_south / u[ south ] ) )
- * ( velocity_x_south * velocity_x_south + velocity_y_south * velocity_y_south ) / 2
- + velocity_y_south * std::sqrt( this->gamma * pressure_south / u[ south ] )
- * ( velocity_y_south + std::sqrt( this->gamma * pressure_south / u[ south ] ) )
- + this->gamma * pressure_south / u[ south ]
- * ( velocity_y_south + std::sqrt( this->gamma * pressure_south / u[ south ] ) )
- / ( this->gamma - 1 )
- )
- - ( u[ center ] / ( 2 * this->gamma ) )
- * (
- ( velocity_y_center - std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- * ( velocity_x_center * velocity_x_center + velocity_y_center * velocity_y_center ) / 2
- + velocity_y_center * std::sqrt( this->gamma * pressure_center / u[ center ] )
- * ( std::sqrt( this->gamma * pressure_center / u[ center ] ) - velocity_y_center )
- + this->gamma * pressure_center / u[ center ]
- * ( velocity_y_center - std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- / ( this->gamma - 1 )
- )
- + ( u[ north ] / ( 2 * this->gamma ) )
- * (
- ( velocity_y_north - std::sqrt( this->gamma * pressure_north / u[ north ] ) )
- * ( velocity_x_north * velocity_x_north + velocity_y_north * velocity_y_north ) / 2
- + velocity_y_north * std::sqrt( this->gamma * pressure_north / u[ north ] )
- * ( std::sqrt( this->gamma * pressure_north / u[ north ] ) - velocity_y_north )
- + this->gamma * pressure_north / u[ north ]
- * ( velocity_y_north - std::sqrt( this->gamma * pressure_north / u[ north ] ) )
- / ( this->gamma - 1 )
- )
+ -hyInverse * (
+ this->positiveEnergyFlux( u[ center ], velocity_y_center, velocity_x_center, 0, pressure_center)
+ - this->positiveEnergyFlux( u[ south ], velocity_y_south , velocity_x_south , 0, pressure_south )
+ - this->negativeEnergyFlux( u[ center ], velocity_y_center, velocity_x_center, 0, pressure_center)
+ + this->negativeEnergyFlux( u[ north ], velocity_y_north , velocity_x_north , 0, pressure_north )
);
}
diff --git a/examples/inviscid-flow-sw/UpwindMomentumBase.h b/examples/inviscid-flow-sw/UpwindMomentumBase.h
index 29d4612c1df34e60f9bf08771673cbde2930c474..2cb3869d074b8e3eafd53bbf186bb647e6882b6a 100644
--- a/examples/inviscid-flow-sw/UpwindMomentumBase.h
+++ b/examples/inviscid-flow-sw/UpwindMomentumBase.h
@@ -87,9 +87,9 @@ class UpwindMomentumBase
if ( machNumber <= -1.0 )
return 0;
else if ( machNumber <= 0.0 )
- return density * speedOfSound * speedOfSound / ( 2 * this->gamma ) * ( machNumber + 1.0 ) * velocity_other / speedOfSound;
+ return density * speedOfSound / ( 2 * this->gamma ) * ( machNumber + 1.0 ) * velocity_other;
else if ( machNumber <= 1.0 )
- return density * speedOfSound * speedOfSound / ( 2 * this->gamma ) * ( ( 2.0 * this->gamma - 1.0 ) * machNumber + 1.0 ) * velocity_other / speedOfSound;
+ return density * speedOfSound / ( 2 * this->gamma ) * ( ( 2.0 * this->gamma - 1.0 ) * machNumber + 1.0 ) * velocity_other;
else
return density * velocity_main * velocity_other;
};
diff --git a/examples/inviscid-flow-sw/UpwindMomentumX.h b/examples/inviscid-flow-sw/UpwindMomentumX.h
index c9cafec4413c6484e41d0968d77bdda463b00c5b..3a6b1abecf3adf5378be5a44d87d4899b6d37d73 100644
--- a/examples/inviscid-flow-sw/UpwindMomentumX.h
+++ b/examples/inviscid-flow-sw/UpwindMomentumX.h
@@ -172,36 +172,16 @@ class UpwindMomentumX< Meshes::Grid< 2, MeshReal, Device, MeshIndex >, Real, Ind
const RealType& velocity_y_south = this->velocity.template getData< DeviceType >()[ 1 ].template getData< DeviceType >()[ south ];
return -hxInverse * (
- u[ center ] / ( 2 * this->gamma )
- * ( 2 * ( this->gamma - 1 ) * velocity_x_center * velocity_x_center
- + ( velocity_x_center + std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- * ( velocity_x_center + std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- )
- - u[ west ] / ( 2 * this->gamma )
- * ( 2 * ( this->gamma - 1 ) * velocity_x_west * velocity_x_west
- + ( velocity_x_west + std::sqrt( this->gamma * pressure_west / u[ west ] ) )
- * ( velocity_x_west + std::sqrt( this->gamma * pressure_west / u[ west ] ) )
- )
- - u[ center ] / ( 2 * this->gamma )
- * ( velocity_x_center - std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- * ( velocity_x_center - std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- + u[ east ] / ( 2 * this->gamma )
- * ( velocity_x_east - std::sqrt( this->gamma * pressure_east / u[ east ] ) )
- * ( velocity_x_east - std::sqrt( this->gamma * pressure_east / u[ east ] ) )
+ this->positiveMainMomentumFlux( u[ center ], velocity_x_center, pressure_center )
+ - this->positiveMainMomentumFlux( u[ west ], velocity_x_west , pressure_west )
+ - this->negativeMainMomentumFlux( u[ center ], velocity_x_center, pressure_center )
+ + this->negativeMainMomentumFlux( u[ east ], velocity_x_east , pressure_east )
)
-hyInverse * (
- u[ center ] / ( 2 * this->gamma )
- * ( ( 2 * this->gamma - 1 ) * velocity_x_center + std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- * velocity_y_center
- - u[ south ] / ( 2 * this->gamma )
- * ( ( 2 * this->gamma - 1 ) * velocity_x_south + std::sqrt( this->gamma * pressure_south / u[ south ] ) )
- * velocity_y_south
- - u[ center ] / ( 2 * this->gamma )
- * ( velocity_x_center - std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- * velocity_y_center
- + u[ north ] / ( 2 * this->gamma )
- * ( velocity_x_north - std::sqrt( this->gamma * pressure_north / u[ north ] ) )
- * velocity_y_north
+ this->positiveOtherMomentumFlux( u[ center ], velocity_x_center, velocity_y_center, pressure_center )
+ - this->positiveOtherMomentumFlux( u[ south ], velocity_x_south , velocity_y_south , pressure_south )
+ - this->negativeOtherMomentumFlux( u[ center ], velocity_x_center, velocity_y_center, pressure_center )
+ + this->negativeOtherMomentumFlux( u[ north ], velocity_x_north , velocity_y_north , pressure_north )
);
diff --git a/examples/inviscid-flow-sw/UpwindMomentumX.h~ b/examples/inviscid-flow-sw/UpwindMomentumX.h~
index 7f72b31daefabb61fe570b56fc7c41b02cbb5f22..095ac28c2160b3c9f26b97c8114ccfefbce8b316 100644
--- a/examples/inviscid-flow-sw/UpwindMomentumX.h~
+++ b/examples/inviscid-flow-sw/UpwindMomentumX.h~
@@ -80,8 +80,8 @@ class UpwindMomentumX< Meshes::Grid< 1, MeshReal, Device, MeshIndex >, Real, Ind
return -hxInverse * (
this->positiveMainMomentumFlux( u[ center ], velocity_x_center, pressure_center )
- this->positiveMainMomentumFlux( u[ west ], velocity_x_west , pressure_west )
- - this->positiveMainMomentumFlux( u[ center ], velocity_x_center, pressure_center )
- + this->positiveMainMomentumFlux( u[ east ], velocity_x_east , pressure_east )
+ - this->negativeMainMomentumFlux( u[ center ], velocity_x_center, pressure_center )
+ + this->negativeMainMomentumFlux( u[ east ], velocity_x_east , pressure_east )
);
}
@@ -172,37 +172,17 @@ class UpwindMomentumX< Meshes::Grid< 2, MeshReal, Device, MeshIndex >, Real, Ind
const RealType& velocity_y_south = this->velocity.template getData< DeviceType >()[ 1 ].template getData< DeviceType >()[ south ];
return -hxInverse * (
- u[ center ] / ( 2 * this->gamma )
- * ( 2 * ( this->gamma - 1 ) * velocity_x_center * velocity_x_center
- + ( velocity_x_center + std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- * ( velocity_x_center + std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- )
- - u[ west ] / ( 2 * this->gamma )
- * ( 2 * ( this->gamma - 1 ) * velocity_x_west * velocity_x_west
- + ( velocity_x_west + std::sqrt( this->gamma * pressure_west / u[ west ] ) )
- * ( velocity_x_west + std::sqrt( this->gamma * pressure_west / u[ west ] ) )
- )
- - u[ center ] / ( 2 * this->gamma )
- * ( velocity_x_center - std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- * ( velocity_x_center - std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- + u[ east ] / ( 2 * this->gamma )
- * ( velocity_x_east - std::sqrt( this->gamma * pressure_east / u[ east ] ) )
- * ( velocity_x_east - std::sqrt( this->gamma * pressure_east / u[ east ] ) )
+ this->positiveMainMomentumFlux( u[ center ], velocity_x_center, pressure_center )
+ - this->positiveMainMomentumFlux( u[ west ], velocity_x_west , pressure_west )
+ - this->negativeMainMomentumFlux( u[ center ], velocity_x_center, pressure_center )
+ + this->negativeMainMomentumFlux( u[ east ], velocity_x_east , pressure_east )
)
- -hyInverse * (
- u[ center ] / ( 2 * this->gamma )
- * ( ( 2 * this->gamma - 1 ) * velocity_x_center + std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- * velocity_y_center
- - u[ south ] / ( 2 * this->gamma )
- * ( ( 2 * this->gamma - 1 ) * velocity_x_south + std::sqrt( this->gamma * pressure_south / u[ south ] ) )
- * velocity_y_south
- - u[ center ] / ( 2 * this->gamma )
- * ( velocity_x_center - std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- * velocity_y_center
- + u[ north ] / ( 2 * this->gamma )
- * ( velocity_x_north - std::sqrt( this->gamma * pressure_north / u[ north ] ) )
- * velocity_y_north
- );
+/* -hyInverse * (
+ this->positiveOtherMomentumFlux( u[ center ], velocity_x_center, velocity_y_center, pressure_center )
+ - this->positiveOtherMomentumFlux( u[ south ], velocity_x_south , velocity_y_south , pressure_south )
+ - this->negativeOtherMomentumFlux( u[ center ], velocity_x_center, velocity_y_center, pressure_center )
+ + this->negativeOtherMomentumFlux( u[ north ], velocity_x_north , velocity_y_north , pressure_north )
+ )*/;
}
diff --git a/examples/inviscid-flow-sw/UpwindMomentumY.h b/examples/inviscid-flow-sw/UpwindMomentumY.h
index 5d3a015c6d0e1927acddc62692e4ce05233d3e45..9859f0618d92d723aba0a4ce3735b00c39e808e5 100644
--- a/examples/inviscid-flow-sw/UpwindMomentumY.h
+++ b/examples/inviscid-flow-sw/UpwindMomentumY.h
@@ -156,36 +156,16 @@ class UpwindMomentumY< Meshes::Grid< 2, MeshReal, Device, MeshIndex >, Real, Ind
const RealType& velocity_y_south = this->velocity.template getData< DeviceType >()[ 1 ].template getData< DeviceType >()[ south ];
return -hxInverse * (
- u[ center ] / ( 2 * this->gamma )
- * ( ( 2 * this->gamma - 1 ) * velocity_y_center + std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- * velocity_x_center
- - u[ west ] / ( 2 * this->gamma )
- * ( ( 2 * this->gamma - 1 ) * velocity_y_west + std::sqrt( this->gamma * pressure_west / u[ west ] ) )
- * velocity_x_west
- - u[ center ] / ( 2 * this->gamma )
- * ( velocity_y_center - std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- * velocity_x_center
- + u[ east ] / ( 2 * this->gamma )
- * ( velocity_y_east - std::sqrt( this->gamma * pressure_east / u[ east ] ) )
- * velocity_x_east
+ this->positiveOtherMomentumFlux( u[ center ], velocity_y_center, velocity_x_center, pressure_center )
+ - this->positiveOtherMomentumFlux( u[ west ], velocity_y_west , velocity_x_west , pressure_west )
+ - this->negativeOtherMomentumFlux( u[ center ], velocity_y_center, velocity_x_center, pressure_center )
+ + this->negativeOtherMomentumFlux( u[ east ], velocity_y_east , velocity_x_east , pressure_east )
)
-hyInverse * (
- u[ center ] / ( 2 * this->gamma )
- * ( 2 * ( this->gamma - 1 ) * velocity_y_center * velocity_y_center
- + ( velocity_y_center + std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- * ( velocity_y_center + std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- )
- - u[ south ] / ( 2 * this->gamma )
- * ( 2 * ( this->gamma - 1 ) * velocity_y_south * velocity_y_south
- + ( velocity_y_south + std::sqrt( this->gamma * pressure_south / u[ south ] ) )
- * ( velocity_y_south + std::sqrt( this->gamma * pressure_south / u[ south ] ) )
- )
- - u[ center ] / ( 2 * this->gamma )
- * ( velocity_y_center - std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- * ( velocity_y_center - std::sqrt( this->gamma * pressure_center / u[ center ] ) )
- + u[ north ] / ( 2 * this->gamma )
- * ( velocity_y_north - std::sqrt( this->gamma * pressure_north / u[ north ] ) )
- * ( velocity_y_north - std::sqrt( this->gamma * pressure_north / u[ north ] ) )
+ this->positiveMainMomentumFlux( u[ center ], velocity_y_center, pressure_center )
+ - this->positiveMainMomentumFlux( u[ south ], velocity_y_south , pressure_south )
+ - this->negativeMainMomentumFlux( u[ center ], velocity_y_center, pressure_center )
+ + this->negativeMainMomentumFlux( u[ north ], velocity_y_north , pressure_north )
);
}