Commit 355716a4 authored by Steffen Schotthöfer's avatar Steffen Schotthöfer
Browse files

refined volume output check; implemented analytical solution for sigma=1 - not yet tested

parent 4c0a7d20
Pipeline #111467 failed with stage
in 15 minutes and 46 seconds
......@@ -10,7 +10,7 @@
% Output directory
OUTPUT_DIR = ../result
% Output file
OUTPUT_FILE = exampleMN_dualM_Analytic
OUTPUT_FILE = exampleMN
% Log directory
LOG_DIR = ../result/logs
% Mesh File
......@@ -28,7 +28,7 @@ CFL_NUMBER = 0.7
% Final time for simulation
TIME_FINAL = 0.3
% Maximal Moment degree
MAX_MOMENT_SOLVER = 1
MAX_MOMENT_SOLVER = 5
%
%% Entropy settings
ENTROPY_FUNCTIONAL = MAXWELL_BOLZMANN
......@@ -60,4 +60,4 @@ QUAD_ORDER = 12
%
% ----- Output ----
%
VOLUME_OUTPUT = (ANALYTIC, MOMENTS, DUAL_MOMENTS)
VOLUME_OUTPUT = ( MINIMAL, MOMENTS, DUAL_MOMENTS)
......@@ -14,8 +14,8 @@ OUTPUT_FILE = examplePN
% Log directory
LOG_DIR = ../result/logs
% Mesh File
%MESH_FILE = linesource.su2
MESH_FILE = linesource_debug.su2
MESH_FILE = linesource.su2
%MESH_FILE = linesource_debug.su2
%
PROBLEM = LINESOURCE
......@@ -25,11 +25,11 @@ PROBLEM = LINESOURCE
% Solver type
SOLVER = PN_SOLVER
% CFL number
CFL_NUMBER = 0.6
CFL_NUMBER = 0.7
% Final time for simulation
TIME_FINAL = 0.3
% Maximal Moment degree
MAX_MOMENT_SOLVER = 2
MAX_MOMENT_SOLVER = 3
% ---- Boundary Conditions ----
% Example: BC_DIRICLET = (dummyMarker1, dummyMarker2)
......
......@@ -20,7 +20,7 @@ PROBLEM = LINESOURCE
% Solver type
SOLVER = SN_SOLVER
% CFL number
CFL_NUMBER = 0.4
CFL_NUMBER = 0.7
% Final time for simulation
TIME_FINAL = 0.3
......@@ -31,6 +31,6 @@ BC_DIRICHLET = ( void )
%
%% Quadrature Specifications
% Quadrature Type
QUAD_TYPE = LEBEDEV
QUAD_TYPE = GAUSS_LEGENDRE_TENSORIZED
% Quadrature Order
QUAD_ORDER = 5
QUAD_ORDER = 6
......@@ -429,10 +429,13 @@ void Config::SetPostprocessing() {
}
// Check, if the choice of volume output is compatible to the solver
std::vector<VOLUME_OUTPUT> supportedGroups;
for( unsigned short idx_volOutput = 0; idx_volOutput < _nVolumeOutput; idx_volOutput++ ) {
switch( _solverName ) {
case SN_SOLVER:
if( _volumeOutput[idx_volOutput] != MINIMAL && _volumeOutput[idx_volOutput] != ANALYTIC ) {
supportedGroups = { MINIMAL, ANALYTIC };
if( supportedGroups.end() == std::find( supportedGroups.begin(), supportedGroups.end(), _volumeOutput[idx_volOutput] ) ) {
ErrorMessages::Error( "SN_SOLVER only supports volume output MINIMAL and ANALYTIC.\nPlease check your .cfg file.",
CURRENT_FUNCTION );
}
......@@ -443,30 +446,34 @@ void Config::SetPostprocessing() {
}
break;
case MN_SOLVER:
if( _volumeOutput[idx_volOutput] != MINIMAL && _volumeOutput[idx_volOutput] != MOMENTS &&
_volumeOutput[idx_volOutput] != DUAL_MOMENTS && _volumeOutput[idx_volOutput] != ANALYTIC ) {
supportedGroups = { MINIMAL, MOMENTS, DUAL_MOMENTS, ANALYTIC };
if( supportedGroups.end() == std::find( supportedGroups.begin(), supportedGroups.end(), _volumeOutput[idx_volOutput] ) ) {
ErrorMessages::Error(
"MN_SOLVER only supports volume output ANALYTIC, MINIMAL, MOMENTS and DUAL_MOMENTS.\nPlease check your .cfg file.",
CURRENT_FUNCTION );
}
if( _volumeOutput[idx_volOutput] == ANALYTIC && _problemName != PROBLEM_LineSource ) {
ErrorMessages::Error(
"Analytical solution (VOLUME_OUTPUT=ANALYTIC) is only available for the PROBLEM=LINESOURCE.\nPlease check your .cfg file.",
CURRENT_FUNCTION );
ErrorMessages::Error( "Analytical solution (VOLUME_OUTPUT=ANALYTIC) is only available for the PROBLEM=LINESOURCE.\nPlease "
"check your .cfg file.",
CURRENT_FUNCTION );
}
break;
case PN_SOLVER:
if( _volumeOutput[idx_volOutput] != MINIMAL && _volumeOutput[idx_volOutput] != MOMENTS && _volumeOutput[idx_volOutput] != ANALYTIC ) {
supportedGroups = { MINIMAL, MOMENTS, ANALYTIC };
if( supportedGroups.end() == std::find( supportedGroups.begin(), supportedGroups.end(), _volumeOutput[idx_volOutput] ) ) {
ErrorMessages::Error( "PN_SOLVER only supports volume output ANALYTIC, MINIMAL and MOMENTS.\nPlease check your .cfg file.",
CURRENT_FUNCTION );
}
if( _volumeOutput[idx_volOutput] == ANALYTIC && _problemName != PROBLEM_LineSource ) {
ErrorMessages::Error(
"Analytical solution (VOLUME_OUTPUT=ANALYTIC) is only available for the PROBLEM=LINESOURCE.\nPlease check your .cfg file.",
CURRENT_FUNCTION );
ErrorMessages::Error( "Analytical solution (VOLUME_OUTPUT=ANALYTIC) is only available for the PROBLEM=LINESOURCE.\nPlease "
"check your .cfg file.",
CURRENT_FUNCTION );
}
break;
case CSD_SN_SOLVER:
supportedGroups = { MINIMAL };
if( _volumeOutput[idx_volOutput] != MINIMAL ) {
ErrorMessages::Error( "CSD_SN_SOLVER only supports volume output MINIMAL.\nPlease check your .cfg file.", CURRENT_FUNCTION );
}
......
......@@ -49,14 +49,12 @@ double LineSource::HelperIntRho_ptc( double R, double t ) {
double LineSource::HelperRho_ptc( double R, double t ) {
double result = HelperRho_ptc1( R, t ) + HelperRho_ptc2( R, t );
return result;
}
double LineSource::HelperRho_ptc1( double R, double t ) {
double gamma = R / t;
double result = exp( -t ) / ( 4.0 * M_PI * R * t ) * log( ( 1 + gamma ) / ( 1 - gamma ) );
return result;
}
......@@ -78,9 +76,10 @@ double LineSource::HelperIntRho_ptc2( double t, double gamma ) {
double u = 0;
double integral = 0;
double beta = 0;
double q = 0;
complex<double> complexPart;
std::complex<double> beta( 0, 0 );
double q = 0;
std::complex<double> complexPart( 0, 0 );
std::complex<double> com_one( 0, 1 );
// compute the integral using the midpoint rule
// integral is from 0 to pi
......@@ -90,17 +89,17 @@ double LineSource::HelperIntRho_ptc2( double t, double gamma ) {
double stepsize = M_PI / (double)numsteps;
u = 0.5 * stepsize;
q = ( 1 + gamma ) / ( 1 - gamma );
q = ( 1.0 + gamma ) / ( 1.0 - gamma );
for( int i = 0; i < numsteps; i++ ) {
u = u + i * stepsize;
u = u + (double)i * stepsize;
// function evaluation
beta = ( log( q ) + _Complex_I * u ) / ( gamma + _Complex_I * tan( 0.5 * u ) );
complexPart = ( gamma + _Complex_I * tan( 0.5 * u ) ) * beta * beta * beta * exp( 0.5 * t * ( 1 - gamma * gamma ) * beta );
beta = ( log( q ) + com_one * u ) / ( gamma + com_one * tan( 0.5 * u ) );
complexPart = ( gamma + com_one * tan( 0.5 * u ) ) * beta * beta * beta * exp( 0.5 * t * ( 1 - gamma * gamma ) * beta );
integral += ( 1 / cos( 0.5 * u ) ) * ( 1 / cos( 0.5 * u ) ) * complexPart.real();
}
return integral;
}
......
......@@ -290,7 +290,7 @@ double MNSolver::WriteOutputFields( unsigned idx_pseudoTime ) {
for( unsigned idx_cell = 0; idx_cell < _nCells; ++idx_cell ) {
double time = idx_pseudoTime * _dE;
double sigma = 0;
double sigma = 1;
_outputFields[idx_group][0][idx_cell] = _problem->GetAnalyticalSolution(
_mesh->GetCellMidPoints()[idx_cell][0], _mesh->GetCellMidPoints()[idx_cell][1], time, sigma );
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment