screen output is now correct mass for pn and mn

Former-commit-id: f449aac2

code/src/common/io.cpp

 /*!
  * \file io.cpp
  * \brief Set of utility io functions for rtns
- * \author S. Schotthoefer, J. Wolters et al.
+ * \author  J. Kusch, S. Schotthoefer, P. Stammer,  J. Wolters, T. Xiao
  */
 
 #include "common/io.h"

code/src/problems/linesource.cpp

@@ -193,7 +193,7 @@ VectorVector LineSource_PN::SetupIC() {
     VectorVector psi( _mesh->GetNumCells(), Vector( ntotalEquations, 0 ) );    // zero could lead to problems?
     VectorVector cellMids = _mesh->GetCellMidPoints();
 
-    // Initial condition is dirac impulse at (x,y) = (0,0) ==> constant in angle ==> all moments are zero.
+    // Initial condition is dirac impulse at (x,y) = (0,0) ==> constant in angle ==> all moments - exept first - are zero.
     double t = 3.2e-4;    // pseudo time for gaussian smoothing (Approx to dirac impulse)
     for( unsigned j = 0; j < cellMids.size(); ++j ) {
         double x  = cellMids[j][0];

code/src/solvers/mnsolver.cpp

@@ -263,6 +263,8 @@ double MNSolver::WriteOutputFields( unsigned idx_pseudoTime ) {
     for( unsigned idx_cell = 0; idx_cell < _nCells; ++idx_cell ) {
         mass += _sol[idx_cell][0] * _areas[idx_cell];    // Should probably go to postprocessing
     }
+    mass *= firstMomentScaleFactor;
+
     if( _settings->GetOutputFrequency() != 0 && idx_pseudoTime % (unsigned)_settings->GetOutputFrequency() == 0 ||
         idx_pseudoTime == _nEnergies - 1 /* need sol at last iteration */ ) {
 

code/src/solvers/pnsolver.cpp

@@ -380,19 +380,21 @@ double PNSolver::WriteOutputFields( unsigned idx_pseudoTime ) {
     for( unsigned idx_cell = 0; idx_cell < _nCells; ++idx_cell ) {
         mass += _sol[idx_cell][0] * _areas[idx_cell];    // Should probably go to postprocessing
     }
+    mass *= firstMomentScaleFactor;
+
     if( _settings->GetOutputFrequency() != 0 && idx_pseudoTime % (unsigned)_settings->GetOutputFrequency() == 0 ||
         idx_pseudoTime == _nEnergies - 1 /* need sol at last iteration */ ) {
         for( unsigned idx_group = 0; idx_group < nGroups; idx_group++ ) {
             switch( _settings->GetVolumeOutput()[idx_group] ) {
                 case MINIMAL:
                     for( unsigned idx_cell = 0; idx_cell < _nCells; ++idx_cell ) {
-                        _outputFields[idx_group][0][idx_cell] = _sol[idx_cell][0];
+                        _outputFields[idx_group][0][idx_cell] = firstMomentScaleFactor * _sol[idx_cell][0];
                     }
                     break;
                 case MOMENTS:
                     for( unsigned idx_sys = 0; idx_sys < _nTotalEntries; idx_sys++ ) {
                         for( unsigned idx_cell = 0; idx_cell < _nCells; ++idx_cell ) {
-                            _outputFields[idx_group][idx_sys][idx_cell] = firstMomentScaleFactor * _sol[idx_cell][idx_sys];
+                            _outputFields[idx_group][idx_sys][idx_cell] = _sol[idx_cell][idx_sys];
                         }
                     }
                     break;