//=============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
//
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE.  See the above copyright notice for more information.
//=============================================================================
#include "smtk/session/rgg/operators/ExportInp.h"

#include "smtk/session/rgg/Resource.h"
#include "smtk/session/rgg/Session.h"

#include "smtk/io/Logger.h"

#include "smtk/PublicPointerDefs.h"
#include "smtk/attribute/Attribute.h"
#include "smtk/attribute/DirectoryItem.h"
#include "smtk/attribute/DoubleItem.h"
#include "smtk/attribute/GroupItem.h"
#include "smtk/attribute/IntItem.h"
#include "smtk/attribute/ModelEntityItem.h"
#include "smtk/attribute/StringItem.h"
#include "smtk/attribute/VoidItem.h"

#include "smtk/io/Logger.h"

#include "smtk/model/AuxiliaryGeometry.h"
#include "smtk/model/Model.h"
#include "smtk/model/Resource.h"
#include "smtk/model/Resource.txx"

#include "smtk/session/rgg/operators/CreateModel.h"

#include "smtk/session/rgg/json/jsonAssembly.h"
#include "smtk/session/rgg/json/jsonCore.h"
#include "smtk/session/rgg/json/jsonDuct.h"
#include "smtk/session/rgg/json/jsonPin.h"

#include "smtk/session/rgg/Assembly.h"
#include "smtk/session/rgg/Core.h"
#include "smtk/session/rgg/Duct.h"
#include "smtk/session/rgg/Pin.h"

#include <boost/cstdint.hpp>
//force to use filesystem version 3
#define BOOST_FILESYSTEM_VERSION 3
#include <boost/filesystem.hpp>

#include <fstream> // For file I/O

#include "smtk/session/rgg/operators/ExportInp_xml.h"

using namespace smtk::model;
using namespace smtk::session::rgg;
using json = nlohmann::json;

namespace
{
void  writeHeader(std::ofstream& output, std::string type)
{
  output << "!   ########################################################\n";
  output << "!   " << type << " File Generated by SMTK RGG Session\n";
  output << "!   ########################################################\n";
}

// Cubit uses material definitions to partition like mesh elements into blocks.
// We would rather handle the material assignment (and mesh attribution)
// ourselves, so we compose "materials" according to the model component.
void writeMaterial(std::ofstream& output, const Assembly& assy, smtk::model::Model model)
{
  // Find out all used materials and preserve their traversal order:
  //    - ducts
  //      - segments
  //        - layers
  //    - pins
  //      - pieces
  //        - layers
  //      - cell material
  std::vector<std::pair<std::string,std::string>> usedMaterials;

  smtk::model::EntityRef ductAux(model.resource(), assy.associatedDuct());
  Duct duct = json::parse(
      ductAux.stringProperty(Duct::propDescription)[0]);
  int segmentIndex = 0;
  const auto& segments = duct.segments();
  for(const auto& segment: segments)
  {
    int layerIndex = 0;
    for (const Duct::Segment::layer& layer : segment.layers)
    {
      if (std::get<0>(layer))
      {
        std::string m = duct.name() + "_" + std::to_string(segmentIndex) +
          "_" + std::to_string(layerIndex);
        usedMaterials.push_back(std::make_pair(m,m));
      }
      ++layerIndex;
    }
    ++segmentIndex;
  }

  const auto& uuidToSchema = assy.layout();
  int pinIndex = 0;
  for (const auto& iter:uuidToSchema)
  {
    smtk::model::EntityRef pinAux(model.resource(), iter.first);
    Pin pin = json::parse(
        pinAux.stringProperty(Pin::propDescription)[0]);
    int pieceIndex = 0;
    const auto& pieces = pin.pieces();
    for (const auto& piece : pieces)
    {
      int layerIndex = 0;
      const auto&  layerMaterials = pin.layerMaterials();
      for (const auto& layerMaterial: layerMaterials)
      {
        if (layerMaterial.subMaterialIndex)
        {
          std::string m = pin.name() + "_" + std::to_string(pinIndex) +
            "_" + std::to_string(pieceIndex) + "_" + std::to_string(layerIndex);
          usedMaterials.push_back(std::make_pair(m,m));
        }
        ++layerIndex;
      }
      ++pieceIndex;
    }
    if (pin.cellMaterialIndex())
    {
      std::string m = pin.name() + "_" + std::to_string(pinIndex);
      usedMaterials.push_back(std::make_pair(m,m));
    }
    ++pinIndex;
  }

  output << "Materials " << usedMaterials.size();
  for (const auto& mAndL : usedMaterials)
  {
    output << " " << mAndL.first << " " << mAndL.second;
  }
  output << "\n";
}

void writeDuct(std::ofstream& output, const Assembly& assy,
                   smtk::model::Model model)
{
  // Format for each segment along z axis - reverse engineered from meshkit source code
  // Duct <number of layers> <X> <Y> <Z0> <Z1> <R0> <R1>... <R0Material> <R1Material>...
  // If it's a rect duct, then its <R0_0> <R0_1>  <R1_0> <R1_1>...
  // <R0Material_0> <R0Material_1> <R1Material_0> <R1Material_1>... after defining z values
  // Ex. a duct which has two segments from 0~5 and 5-10 and thickness as 10.
  // the first layer has 1 radius as 1 and the second has 2 radius as 0.6 and 1.
  // Dimensions 1 0.00000 0.00000 0.00000 5.00000 10.0000 Unknown
  // Dimensions 2 0.00000 0.00000 5.00000 10.0000 6.00000 10.0000 Unknown Unknown
  bool isHex = static_cast<Core::GeomType>(model.integerProperty(Core::geomDescription)[0]) ==
      Core::GeomType::Hex;

  auto coreGroup = model.resource()->findEntitiesByProperty("rggType", Core::typeDescription)[0];
  Core core = json::parse(coreGroup.stringProperty(Core::propDescription)[0]);
  std::pair<double, double> ductThickness = core.ductThickness();

  smtk::model::EntityRef ductAux(model.resource(), assy.associatedDuct());
  Duct duct = json::parse(
      ductAux.stringProperty(Duct::propDescription)[0]);
  int x{0}, y{0};
  const auto& segments = duct.segments();
  int segmentIndex = 0;
  for (const auto& segment : segments)
  {
    output << "Duct " << segment.layers.size() << " ";
    output << x << " " << y << " "; // For now rgg session does not allow custom center coordinates.
    output << segment.baseZ << " " << segment.baseZ + segment.height << " ";
    const auto& layers = segment.layers;
    // Export radius info
    for (const auto& layer : layers)
    {
      if (isHex)
      {
        output << ductThickness.first * std::get<1>(layer) << " ";
      }
      else
      {
        output << ductThickness.first * std::get<1>(layer) << " "
               << ductThickness.second * std::get<2>(layer) << " ";
      }
    }
    // Export material info
    int layerIndex = 0;
    for (const auto& layer : layers)
    {
      if (isHex)
      {
        output << duct.name() << "_" << std::to_string(segmentIndex)
               << "_" << std::to_string(layerIndex) << " ";
      }
      else
      {
        output << duct.name() << "_" << std::to_string(segmentIndex)
               << "_" << std::to_string(layerIndex) << " ";
        output << duct.name() << "_" << std::to_string(segmentIndex)
               << "_" << std::to_string(layerIndex) << " ";
      }
      ++layerIndex;
    }
    ++segmentIndex;
    output << "\n";
  }
}

void writePins(std::ofstream& output, const Assembly& assy,
                   smtk::model::Model model)
{
  // Format for each segment along z axis - reverse engineered from meshkit source code
  // cylinder/frustrum <number of layers> <X> <Y> <Z0> <Z1> <R0> <R1>... <R0Material> <R1Material>...
  // The geometry type of the core does not affect pin
  // Ex. a pin which has three segments from 0~5, 5-10 and 10-20 with material normRadius as 0.5 and 1.
  // Radius for all cylinders are 0.5.
  // PinCell0 PC0 3
  // cylinder 2 0.00000 0.00000 0.00000 5.00000 0.250000 0.500000 activecore coolant
  // cylinder 2 0.00000 0.00000 5.00000 10.0000 0.250000 0.500000 activecore coolant
  // cylinder 2 0.00000 0.00000 10.0000 20.0000 0.250000 0.500000 activecore coolant
  bool isHex = static_cast<Core::GeomType>(model.integerProperty(Core::geomDescription)[0]) ==
      Core::GeomType::Hex;
  std::pair<double, double> pitch = assy.pitch();

  const auto& uuidToSchema = assy.layout();
  output << "pincells " << uuidToSchema.size() << " " << pitch.first;
  if (!isHex)
  {
    output << " " << pitch.second << "\n";
  }
  else
  {
    output << "\n";
  }
  int pinIndex = 0;
  for (const auto& iter:uuidToSchema)
  {
    smtk::model::EntityRef pinAux(model.resource(), iter.first);
    Pin pin = json::parse(
        pinAux.stringProperty(Pin::propDescription)[0]);
    const auto& pieces = pin.pieces();
    const auto&  layerMaterials = pin.layerMaterials();
    int x{0}, y{0};

    output << pin.name() << " " << pin.label() << " " << pieces.size() + (pin.cellMaterialIndex() > 0)<< std::endl;

    double baseZ = pin.zOrigin();
    int pieceIndex = 0;
    for (const auto& piece: pieces)
    {
      bool isCylinder = piece.pieceType == Pin::PieceType::CYLINDER;
      output << (isCylinder ? "cylinder " : "frustum ") << layerMaterials.size() << " ";
      output << x << " " << y << " ";
      output << baseZ << " " << baseZ + piece.length << " ";
      baseZ += piece.length;
      // radius and materials
      for (const auto& layerMaterial: layerMaterials)
      {
        // TODO: Here I just copy the logic from RGG application. Don't know if its right since meshkit example does
        // not include this case
        if (isCylinder)
        {
          output << layerMaterial.normRadius * piece.baseRadius << " ";
        }
        else
        {
          output << layerMaterial.normRadius * piece.baseRadius << " " << layerMaterial.normRadius * piece.topRadius << " ";
        }
      }
      int layerIndex = 0;
      for (const auto& layerMaterial: layerMaterials)
      {
        output << pin.name() << "_" << std::to_string(pinIndex)
               << "_" << std::to_string(pieceIndex) << "_" << std::to_string(layerIndex) << " ";
        ++layerIndex;
      }
      output << "\n";
      ++pieceIndex;
    }
    // FIXME: How to detect if there is no material when user provides custom materials?
    std::string cellM = CreateModel::getMaterial(static_cast<size_t>(pin.cellMaterialIndex()),
                                               model);
    if (cellM != "NoCellMaterial")
    {
      output << "CellMaterial " << pin.zOrigin() << " " << baseZ << " "
             << pin.name() << "_" << std::to_string(pinIndex) << "\n";
    }
    ++pinIndex;
  }
}

// usingAmp: ‘&’ sign can be used for continuing the line input on the next line in meshkit
void writeLattice(std::ofstream& output, const std::pair<int, int>& latticeSize,
                  const Assembly::UuidToSchema uuidToSchema,
                  smtk::model::Model model, std::string keyword, bool useAmp)
{
  bool isHex = static_cast<Core::GeomType>(model.integerProperty(Core::geomDescription)[0]) ==
      Core::GeomType::Hex;
  output << keyword << " ";
  // When its type is rect, y size should be write out first
  if (!isHex)
  {
    output << latticeSize.second << " ";
  }
  output << latticeSize.first << std::endl;;

  // Create a [schemaCoord]=label map
  std::map<std::pair<int,int>, std::string> coordToLabel;
  for (const auto& iter: uuidToSchema)
  {
    smtk::model::EntityRef pinAux(model.resource(), iter.first);
    if (!pinAux.hasStringProperty("label"))
    {
      smtkErrorMacro(smtk::io::Logger::instance(), "encounter a pin does not have label string property");
      continue;
    }
    std::string label = pinAux.stringProperty("label")[0];
    for (const auto& coord: iter.second)
    {
      coordToLabel[coord]=label;
    }
  }

  // The index order of each layer is clockwise, starting from upper left corner of the hex.
  // It's prefined in old RGG application.
  int numRings = latticeSize.first;
  if (isHex && numRings)
  {
    // Initialize the string array
    int numRows = 2*numRings - 1;
    std::vector<std::vector<std::string> > hexArray;
    hexArray.resize(static_cast<size_t>(numRows));
    int numCols = 0;
    int delta=0;
    for(int i = 0; i < numRows; i++)
    {
      if(i<numRings) // top half of the hex lattice(including the middle row)
      {
        numCols = i + numRings;
      }
      else // second half of the hex lattice
      {
        delta++;
        numCols = numRows - delta;
      }
      hexArray[static_cast<size_t>(i)].resize(static_cast<size_t>(numCols));
      std::fill(hexArray[static_cast<size_t>(i)].begin(),
          hexArray[static_cast<size_t>(i)].end(), "XX");
    }

    // Fill in from out ring to inner ring
    for (int k = static_cast<int>(numRings) -1; k >=0; k--)
    {
      int numRowsC{2*k + 1}, startRowC{numRings-1-k};
      int startColC{startRowC}, layerIndexC{0};
      // For each ring, fill it from top row to the lower row
      for (int i = startRowC; i < startRowC + numRowsC;i++)
      { // First or last row
        if (i == startRowC || i == startRowC + numRowsC -1)
        { // Loop through the columns of current row
          for (int j = startColC, offset = 0; j < k+1+startColC; j++, offset++)
          {
            layerIndexC = (i==startRowC) ? offset : 4*k-offset;
            if (coordToLabel.find(std::make_pair(k, layerIndexC)) != coordToLabel.end())
            {
              hexArray[static_cast<size_t>(i)][static_cast<size_t>(j)] =
                                  coordToLabel[std::make_pair(k, layerIndexC)];
            }
          }
        }
        else
        {
          // Fill in the left value
          layerIndexC = 6*k - (i-startRowC);
          if (coordToLabel.find(std::make_pair(k, layerIndexC)) != coordToLabel.end())
          {
            hexArray[static_cast<size_t>(i)][static_cast<size_t>(startColC)] =
                                coordToLabel[std::make_pair(k, layerIndexC)];
          }
          // Fill in the right value
          layerIndexC = k + (i-startRowC);
          int endColC = static_cast<int>(hexArray[i].size()) - 1 - startColC;
          if (coordToLabel.find(std::make_pair(k, layerIndexC)) != coordToLabel.end())
          {
            hexArray[static_cast<size_t>(i)][static_cast<size_t>(endColC)] =
                                coordToLabel[std::make_pair(k, layerIndexC)];
          }
        }
      }
    }
    for (size_t i = 0; i < hexArray.size();i++)
    {
      if (static_cast<int>(i) <= numRings-1)
      {
        size_t numEmptySpace = (static_cast<size_t>(numRings) -1 - i) * 2;
        output << std::string(numEmptySpace, ' ') ;
      }
      else
      {
        size_t numEmptySpace = (i - (static_cast<size_t>(numRings)-1)) * 2;
        output << std::string(numEmptySpace, ' ') ;
      }
      for (size_t j = 0; j < hexArray[i].size();j++)
      {
        output << hexArray[i][j] << "  ";
      }
      if(i < hexArray.size()-1 && useAmp)
      {
        output << "&";
      }
      output << "\n";
    }
  }
  else
  {
    // TODO: add support for rect layout
  }
}

void writeAssemblies(std::ofstream& output, const Core& core,
                   smtk::model::Model model)
{
  bool isHex = static_cast<Core::GeomType>(model.integerProperty(Core::geomDescription)[0]) ==
      Core::GeomType::Hex;
  const auto& layout = core.layout();
  output << "assemblies " << layout.size() << " " <<  core.ductThickness().first << " ";
  if (isHex)
  {
    output << std::endl;
  }
  else
  {
    output << core.ductThickness().second << std::endl;
  }
  for (const auto& iter : layout)
  {
    smtk::model::EntityRef assy(model.resource(), iter.first);
    if (!assy.isValid() || !assy.hasStringProperty("label"))
    {
      continue;
    }
    output << assy.name() << ".exo " << assy.stringProperty("label")[0] << std::endl;
  }
}

}

namespace smtk
{
namespace session
{
namespace rgg
{

ExportInp::Result ExportInp::operateInternal()
{
  EntityRefArray entities = this->parameters()->associatedModelEntities<EntityRefArray>();
  if (entities.empty() || (!entities[0].isAuxiliaryGeometry() && !entities[0].isModel()))
  {
    smtkErrorMacro(this->log(), "Cannot edit a non auxiliary geometry nor model");
    return this->createResult(smtk::operation::Operation::Outcome::FAILED);
  }
  smtk::model::Model model = entities[0].as<smtk::model::Model>();
  smtk::model::ResourcePtr resource = model.resource();

  // Get the export dir
  smtk::attribute::DirectoryItemPtr pathItem = this->parameters()->findDirectory("directory");
  if (pathItem != nullptr && !pathItem->value(0).empty())
  {
    m_dir = pathItem->value(0);
  }
  // TODO: Add verification logic to the path

  // Fetch rgg core and its assemblies
  auto coreGroup = resource->findEntitiesByProperty("rggType", Core::typeDescription)[0];
  auto assemblyGroups = resource->findEntitiesByProperty("rggType", Assembly::typeDescription);


  // Write the commmon.inp file
  if (assemblyGroups.size())
  {
    this->ExportCommonInp(assemblyGroups[0]);
  }
  else
  {
    this->ExportCommonInp(smtk::model::EntityRef());
  }

  // Export each assembly to its inp file
  std::for_each(assemblyGroups.begin(), assemblyGroups.end(), [this](const smtk::model::EntityRef& assy)
  {this->ExportAssembly(assy);});

  // Export core
  this->ExportCore(coreGroup);

  Result result = this->createResult(smtk::operation::Operation::Outcome::SUCCEEDED);

  return result;
}


const char* ExportInp::xmlDescription() const
{
  return ExportInp_xml;
}

void ExportInp::ExportCommonInp(smtk::model::EntityRef assyGroup)
{
  std::string path = m_dir + "/" + "common.inp";
  std::ofstream file(path.c_str());
  if (!file.is_open())
  {
    smtkWarningMacro(smtk::io::Logger::instance(),
                     "Fail to open common.inp file. Abort writing common.inp file");
    return;
  }
  if (assyGroup.isValid())
  {
    Assembly assy = json::parse(
        assyGroup.stringProperty(Assembly::propDescription)[0]);
    const auto& aep = assy.exportParams();
    if (!aep.GeomEngine.empty())
    {
      file << "GeomEngine " << aep.GeomEngine << std::endl;
    }
    if (aep.StartPinId >= 0)
    {
      file << "StartPinId " << aep.StartPinId << std::endl;
    }
    if (!aep.MeshType.empty())
    {
      file << "MeshType " << aep.MeshType << std::endl;
    }
    if (!aep.Info.empty())
    {
      file << "Info " << aep.Info << std::endl;
    }
    if (std::get<0>(aep.HBlock) >= 0)
    {
      file << "HBlock " << std::get<0>(aep.HBlock) << " "
           << std::get<1>(aep.HBlock) << " "
           << std::get<2>(aep.HBlock) << std::endl;
    }
    if (!aep.GeometryType.empty())
    {
      file << "GeometryType " << aep.GeometryType << std::endl;
    }
    if (!aep.Geometry.empty())
    {
      file << "Geometry " << aep.Geometry << std::endl;
    }
    if (!aep.CreateSideset.empty())
    {
      file << "CreateSideset " << aep.CreateSideset << std::endl;
    }
    if (!aep.CreateFiles.empty())
    {
      file << "CreateFiles " << aep.CreateFiles << std::endl;
    }
    if (!aep.SaveExodus.empty())
    {
      file << "SaveExodus " << aep.SaveExodus << std::endl;
    }
    if (aep.MergeTolerance >= 0)
    {
      file << "MergeTolerance " << aep.MergeTolerance << std::endl;
    }
    if (aep.RadialMeshSize >= 0)
    {
      file << "RadialMeshSize " << aep.RadialMeshSize << std::endl;
    }
    if (aep.TetMeshSize >= 0)
    {
      file << "TetMeshSize " << aep.TetMeshSize << std::endl;
    }
    if (aep.AxialMeshSize >= 0)
    {
      file << "AxialMeshSize " << aep.AxialMeshSize << std::endl;
    }
    if (aep.EdgeInterval >= 0)
    {
      file << "EdgeInterval " << aep.EdgeInterval << std::endl;
    }
    if (!aep.MeshScheme.empty())
    {
      file << "MeshScheme " << aep.MeshScheme <<std::endl;
    }
  }
  file << "END\n";
  file.close();
}

void ExportInp::ExportAssembly(smtk::model::EntityRef assyGroup)
{
  if (!assyGroup.isValid())
  {
    smtkErrorMacro(smtk::io::Logger::instance(),
                   "Cannot export an invalid assembly to inp format");
    return;
  }
  Assembly assy = json::parse(
      assyGroup.stringProperty(Assembly::propDescription)[0]);
  const auto& aep = assy.exportParams();

  std::string path = m_dir +'/' + assy.name()+".inp";
  std::ofstream file(path.c_str());
  if(!file.is_open())
  {
    return;
  }
  writeHeader(file, "Assembly");
  // Materials
  file << "!Materials info" <<std::endl;
  writeMaterial(file, assy, assyGroup.owningModel());
  // Duct
  file << "!Duct info" <<std::endl;
  writeDuct(file, assy, assyGroup.owningModel());
  // Pins
  file << "!Pins info" <<std::endl;
  writePins(file, assy, assyGroup.owningModel());
  // lattice
  file << "!Lattice info" <<std::endl;
  writeLattice(file, assy.latticeSize(), assy.layout(), assyGroup.owningModel(), "Assembly", false);

  // Other export info
  file << "!Meshing parameters info" <<std::endl;
  if (!aep.Rotate.empty())
  {
    file << "Rotate " << aep.Rotate << std::endl;
  }
  if (aep.Center)
  {
    file << "Center" << std::endl;
  }
  if (!aep.Section.empty())
  {
    file << "Section " << aep.Section << std::endl;
  }
  if (!aep.Move.empty())
  {
    file << "Move " << aep.Move << std::endl;
  }
  if (aep.NeumannSet_StartId >= 0)
  {
    file << "NeumannSet_StartId " << aep.NeumannSet_StartId << std::endl;
  }
  if (aep.MaterialSet_StartId >= 0)
  {
    file << "MaterialSet_StartId " << aep.MaterialSet_StartId << std::endl;
  }
  if (!aep.Superblocks.empty())
  {
    file << "Superblocks " << aep.Superblocks << std::endl;
  }
  if (aep.NumSuperBlocks >= 0)
  {
    file << "NumSuperBlocks " << aep.NumSuperBlocks << std::endl;
  }
  if (aep.List_MaterialSet_StartId >= 0)
  {
    file << "List_MaterialSet_StartId " << aep.List_MaterialSet_StartId << std::endl;
  }
  if (aep.List_NeumannSet_StartId >= 0)
  {
    file << "List_NeumannSet_StartId " << aep.List_NeumannSet_StartId << std::endl;
  }
  file << "save_exodus\n";
  file << "END ! This is a must\n";
  file.close();
}

void ExportInp::ExportCore(smtk::model::EntityRef coreGroup)
{
  if (!coreGroup.isValid())
  {
    smtkErrorMacro(smtk::io::Logger::instance(),
                   "Cannot export an invalid assembly to inp format");
    return;
  }
  Core core = json::parse(
      coreGroup.stringProperty(Core::propDescription)[0]);

  std::string path = m_dir +'/' + core.name()+".inp";
  std::ofstream file(path.c_str());
  if(!file.is_open())
  {
    return;
  }
  writeHeader(file, "Core");
  // Parameters
  const auto& cep = core.exportParams();
  if (!cep.Geometry.empty())
  {
    file << "Geometry " << cep.Geometry << std::endl;
  }
  if (!cep.GeometryType.empty())
  {
    file << "GeometryType " << cep.GeometryType << std::endl;
  }
  if (!cep.GeomEngine.empty())
  {
    file << "GeomEngine " << cep.GeomEngine << std::endl;
  }
  if (!cep.Extrude.empty())
  {
    file << "Extrude " << cep.Extrude << std::endl;
  }
  if (cep.MergeTolerance >= 0)
  {
    file << "MergeTolerance " << cep.MergeTolerance << std::endl;
  }
  if (!cep.ProblemType.empty())
  {
    file << "ProblemType " << cep.ProblemType << std::endl;
  }
  if (!cep.SaveParallel.empty())
  {
    file << "SaveParallel " << cep.SaveParallel << std::endl;
  }
  if (!cep.Info.empty())
  {
    file << "Info " << cep.Info << std::endl;
  }
  if (!cep.MeshINFO.empty())
  {
    file << "MeshINFO " << cep.MeshINFO << std::endl;
  }
  if (!cep.Symmetry.empty())
  {
    file << "Symmetry " << cep.Symmetry << std::endl;
  }

  // Assembly info
  writeAssemblies(file, core, coreGroup.owningModel());
  // lattice
  file << "!Lattice info" <<std::endl;
  writeLattice(file, core.latticeSize(), core.layout(), coreGroup.owningModel(), "Lattice", true);
  // Somehow if we export these parameters before lattice, meshkit will crash
  if (!cep.BackgroundInfo.first.empty() &&
      !cep.BackgroundInfo.second.empty())
  {
    file << "Background " << cep.BackgroundInfo.first << std::endl;
    ::boost::filesystem::copy_file(cep.BackgroundInfo.second, m_dir + "/" + cep.BackgroundInfo.first);
  }
  if (!cep.NeumannSet.empty())
  {
    for (const auto& neu: cep.NeumannSet)
    {
      file << "NeumannSet " << neu << std::endl;
    }
  }
  if (!cep.OutputFileName.empty())
  {
    file << "OutputFileName " << cep.OutputFileName << std::endl;
  }
  file << "END ! This is a must\n";
  file.close();
}

} // namespace rgg
} //namespace session
} // namespace smtk
