/***********************************************************************
*
* ELMER, A Computational Fluid Dynamics Program.
*
* Copyright 1st April 1995 - , Center for Scientific Computing,
* Finland.
*
* All rights reserved. No part of this program may be used,
* reproduced or transmitted in any form or by any means
* without the written permission of CSC.
*
* Address: Center for Scientific Computing
* Tietotie 6, P.O. BOX 405
* 02101 Espoo, Finland
* Tel. +358 0 457 2001
* Telefax: +358 0 457 2302
* EMail: Jari.Jarvinen@csc.fi
************************************************************************/
/***********************************************************************
Program: ELMER Front
Module: ecif_modelOutputManager.cpp
Language: C++
Date: 13.04.00
Version: 1.00
Author(s): Martti Verho
Revisions:
Abstract: Implementation
************************************************************************/
#include <eio_api.h>
#include "ecif_body.h"
#include "ecif_body2D.h"
#include "ecif_body3D.h"
#include "ecif_bodyLayer.h"
#include "ecif_bodyElement.h"
#include "ecif_bodyElement1D.h"
#include "ecif_bodyElement2D.h"
#include "ecif_bodyElement3D.h"
#include "ecif_bodyElementLoop.h"
#include "ecif_bodyElementGroup.h"
#include "ecif_boundaryCondition.h"
#include "ecif_boundbox.h"
#include "ecif_control.h"
#include "ecif_input.h"
#include "ecif_mesh.h"
// #include "ecif_modelOutputManager.h"
#include "ecif_model_aux.h"
#include "ecif_parameter.h"
#include "ecif_parameterField.h"
#include "ecif_timer.h"
#include "ecif_userinterface.h"
Control* ModelOutputManager::theControlCenter = NULL;
Model* ModelOutputManager::model = NULL;
ModelOutputManager::ModelOutputManager()
{
}
ModelOutputManager::~ModelOutputManager()
{
}
//**************************************************************
//*** EMF output function (Elmer Front's model file format) ***
//**************************************************************
// Save Elmer Front Model File (emf-file)
//
ostream&
ModelOutputManager::emf_output(ostream& out, char* filename)
{
const ModelInfo* mi = model->getModelInfo();
// Info section as comments
out << "!ElmerFront model file" << endl;
out << "!Saved = " << mi->modified << endl;
out << "!Case = " << mi->modelName << " " << mi->problemName << endl;
out << "!Model dir = " << mi->modelDirectory_absolute << endl;
out << "!Include path = " << mi->includePath_absolute << endl;
out << "!Results dir = " << mi->resultsDirectory_absolute << endl;
out << "!Log dir = " << mi->temporaryFilesDirectory_absolute << endl;
out << endl;
emf_outputHeader(out);
emf_outputTimestamps(out);
emf_outputStatistics(out);
emf_outputParameters(out, ECIF_MODEL_PARAMETER);
emf_outputParameters(out, ECIF_SIMULATION_PARAMETER);
emf_outputParameters(out, ECIF_SOLVER_CONTROL);
emf_outputParameters(out, ECIF_CONSTANT);
emf_outputParameters(out, ECIF_COORDINATE);
emf_outputParameters(out, ECIF_DATAFILE);
emf_outputParameters(out, ECIF_EQUATION_VARIABLE);
emf_outputParameters(out, ECIF_EQUATION);
emf_outputParameters(out, ECIF_SOLVER);
emf_outputParameters(out, ECIF_CALCULATOR);
emf_outputParameters(out, ECIF_TIMESTEP);
emf_outputVertexTable(out);
emf_outputElements(out);
emf_outputElementGroups(out);
emf_outputElementLoops(out);
emf_outputBodies(out);
emf_outputParameters(out, ECIF_BODY_PARAMETER);
emf_outputParameters(out, ECIF_BOUNDARY_PARAMETER);
emf_outputParameters(out, ECIF_INITIAL_CONDITION);
emf_outputBoundaryConditions(out);
emf_outputParameters(out, ECIF_MATERIAL);
emf_outputParameters(out, ECIF_BODY_FORCE);
emf_outputParameters(out, ECIF_GRID_PARAMETER);
emf_outputParameters(out, ECIF_GRID_H);
emf_outputParameters(out, ECIF_USER_SETTING);
emf_outputEof(out);
UserInterface* gui = theControlCenter->getGui();
if (mi->gebhardtFactorsNeedsUpdate)
gui->setNeedsUpdate("Gebhardt Factors");
if (mi->meshNeedsUpdate)
gui->setNeedsUpdate("Mesh");
if (mi->solverNeedsUpdate)
gui->setNeedsUpdate("Solver");
return out;
}
// Output bodies.
ostream&
ModelOutputManager::emf_outputBodies(ostream& out)
{
int index = 0;
while (true) {
Body* body = model->getBody(index++);
if (body==NULL) break;
emf_outputSectionStart(out);
body->output_emf(out, ESF_INDENT_SIZE, 0);
emf_outputSectionEnd(out);
}
return out;
}
// Output all boundary conditions (outer, inner).
// NOTE Boundary conditions need extra checking, they
// are not output as standard parameters
ostream&
ModelOutputManager::emf_outputBoundaryConditions(ostream& out)
{
int index = 0;
while (true) {
Parameter* cond = model->getParameter(index++, ECIF_BOUNDARY_CONDITION);
if (cond==NULL) break;
// Output
emf_outputSectionStart(out);
cond->output_emf(out, ESF_INDENT_SIZE, 0, true, true, false);
emf_outputSectionEnd(out);
// Reset
cond->resetValue();
}
return out;
}
// Output boundary groups.
ostream&
ModelOutputManager::emf_outputElementGroups(ostream& out)
{
int index = 0;
while (true) {
BodyElementGroup* beg = model->getBodyElementGroup(index++);
if (beg==NULL) break;
if ( IMPLICIT_GROUP == beg->getGroupType() ) continue;
emf_outputSectionStart(out);
beg->output_emf(out, ESF_INDENT_SIZE, 0);
emf_outputSectionEnd(out);
}
return out;
}
// Output all model elements
ostream&
ModelOutputManager::emf_outputElements(ostream& out)
{
emf_outputVertices(out);
emf_outputEdges(out);
emf_outputFaces(out);
return out;
}
// Output all model geometry vertices
ostream&
ModelOutputManager::emf_outputVertices(ostream& out)
{
int index = 0;
while (true) {
BodyElement* v = model->getVertex(index++);
if (v==NULL) break;
emf_outputSectionStart(out);
// Do not output geometry for a vertex given in the model's vertex-table
if ( model->isInVertexTable(v->Id()) ) {
v->output_emf(out, ESF_INDENT_SIZE, 0, false);
} else {
v->output_emf(out, ESF_INDENT_SIZE, 0, true);
}
emf_outputSectionEnd(out);
}
return out;
}
// Output all model geometry edges
ostream&
ModelOutputManager::emf_outputEdges(ostream& out)
{
const ModelInfo* mi = model->getModelInfo();
double start[3] = {0.0, 0.0, 0.0};
double end1[3] = {0.0, 0.0, 0.0};
double end2[3] = {0.0, 0.0, 0.0};
bool checkSymmetry = model->getSymmetryAxis(start, end1, end2) && (mi->dimension == ECIF_2D);
GcPoint start_p(start);
GcPoint end_p1(end1);
GcPoint end_p2(end1);
int index = 0;
while (true) {
BodyElement* e = model->getEdge(index++);
if (e==NULL) break;
//-Check status
beStatus e_stat = e->getStatus();
if ( e_stat & (BE_DEVIDED | BE_SWAPPED) ) {
continue;
}
emf_outputSectionStart(out);
e->output_emf(out, ESF_INDENT_SIZE, 0);
emf_outputSectionEnd(out);
}
return out;
}
// Output all model geometry faces
ostream&
ModelOutputManager::emf_outputFaces(ostream& out)
{
const ModelInfo* mi = model->getModelInfo();
double start[3] = {0.0, 0.0, 0.0};
double end1[3] = {0.0, 0.0, 0.0};
double end2[3] = {0.0, 0.0, 0.0};
bool checkSymmetry = model->getSymmetryAxis(start, end1, end2) && (mi->dimension == ECIF_3D);
GcPoint start_p(start);
GcPoint end_p1(end1);
GcPoint end_p2(end1);
int index = 0;
while (true) {
BodyElement* f = model->getFace(index++);
if (f==NULL) break;
//-Check status
beStatus f_stat = f->getStatus();
if ( f_stat & (BE_DEVIDED | BE_SWAPPED) ) {
continue;
}
emf_outputSectionStart(out);
f->output_emf(out, ESF_INDENT_SIZE, 0);
emf_outputSectionEnd(out);
}
return out;
}
#if 0
// Output all model elements
ostream&
ModelOutputManager::emf_outputElements(ostream& out)
{
double start[3] = {0.0, 0.0, 0.0};
double end1[3] = {0.0, 0.0, 0.0};
double end2[3] = {0.0, 0.0, 0.0};
bool checkSymmetry = model->getSymmetryAxis(start, end1, end2);
GcPoint start_p(start);
GcPoint end_p1(end1);
GcPoint end_p2(end1);
BodyElement* be = model->getFirstBoundary();
if (be == NULL)
return out;
while (be) {
//-Check status
beStatus be_stat = be->getStatus();
if ( be_stat & (BE_DEVIDED | BE_SWAPPED) ) {
be = getNextBoundary();
continue;
}
// We check if outer boundary is aint possible symmetry axis (2D) or
// symmetry plane (3D).
bool isOnSymmetry = false;
if (checkSymmetry) {
if (modelInfo->dimension == ECIF_2D)
isOnSymmetry = be->isOnSameAxis(start_p, end_p1);
else if (modelInfo->dimension == ECIF_3D)
isOnSymmetry = be->isOnSamePlane(start_p, end_p1, end_p2);
}
emf_outputSectionStart(out);
be->output_emf(out, ESF_INDENT_SIZE, 0, isOnSymmetry);
emf_outputSectionEnd(out);
be = model->getNextBoundary();
}
return out;
}
#endif
// Output all model elements loops
ostream&
ModelOutputManager::emf_outputElementLoops(ostream& out)
{
int index = 0;
while (true) {
BodyElementLoop* bel = model->getBodyElementLoop(index++);
if (bel==NULL) break;
emf_outputSectionStart(out);
bel->output_emf(out, ESF_INDENT_SIZE, 0);
emf_outputSectionEnd(out);
}
return out;
}
// Mark End Of File.
ostream&
ModelOutputManager::emf_outputEof(ostream& out)
{
emf_outputSectionStart(out);
out << "!End Of File" << endl;
return out;
}
// Put model header data into ouput file.
ostream&
ModelOutputManager::emf_outputHeader(ostream& out)
{
char* QM = "\""; // quote-mark
char* ES = "\"\""; // empty string
short IS = ESF_INDENT_SIZE;
const ModelInfo* mi = model->getModelInfo();
const ParallelInfo* pi = model->getParallelInfo();
emf_outputSectionStart(out);
// Header
LibFront::output_string(out, IS, 0, EMF_HEADER, true);
LibFront::output_scalar(out, IS, 1, EMF_CREATED, NULL, mi->created, true);
LibFront::output_scalar(out, IS, 1, EMF_MODIFIED, NULL, mi->modified, true);
if ( mi->hasUserDefinitions ) {
LibFront::output_scalar(out, IS, 1, EMF_HAS_USER_DEFINITIONS, NULL, mi->hasUserDefinitions);
}
// Input version number for version control
// ========================================
//--Store input version number if it is different from the current program version
if ( mi->frontInputVersionNbr > 0 &&
mi->frontInputVersionNbr != mi->frontVersionNbr
) {
LibFront::output_scalar(out, IS, 1, EMF_ELMER_FRONT_INPUT_VERSION, NULL, mi->frontInputVersionNbr);
//--Otherwise store previous input version number if it was given in the input
} else if (mi->frontPreviousInputVersionNbr > 0 ) {
LibFront::output_scalar(out, IS, 1, EMF_ELMER_FRONT_INPUT_VERSION, NULL, mi->frontPreviousInputVersionNbr);
}
LibFront::output_scalar(out, IS, 1, EMF_ELMER_FRONT_VERSION, NULL, mi->frontVersionNbr);
LibFront::output_scalar(out, IS, 1, EMF_TIMESTAMP, NULL, mi->modelFileTs, true);
LibFront::output_scalar(out, IS, 1, EMF_MODEL_STATUS, NULL, mi->modelStatus);
LibFront::output_scalar(out, IS, 1, EMF_MODEL_SOURCE_TYPE, NULL, mi->modelSourceType);
bool cad_src = false;
bool msh_src = false;
switch (mi->modelSourceType) {
case ECIF_CAD_FILE:
cad_src = true; break;
case ECIF_MESH_FILE:
msh_src = true; break;
case ECIF_CAD_AND_MESH_FILE:
cad_src = msh_src = true; break;
}
if (cad_src) {
LibFront::output_scalar(out, IS, 1, EMF_CAD_SOURCE_FILE, NULL, mi->cadSourceFile, true);
}
if (msh_src) {
LibFront::output_scalar(out, IS, 1, EMF_MESH_SOURCE_FILE, NULL, mi->meshSourceFile, true);
}
LibFront::output_scalar(out, IS, 1, EMF_MESH_RESULT_FILE, NULL, mi->meshResultFile, true);
LibFront::output_scalar(out, IS, 1, EMF_MODEL_NAME, NULL, mi->modelName, true);
LibFront::output_scalar(out, IS, 1, EMF_PROBLEM_NAME, NULL, mi->problemName, true);
LibFront::output_scalar(out, IS, 1, EMF_MODEL_DESCRIPTION, NULL, mi->modelDescription, true);
LibFront::output_scalar(out, IS, 1, EMF_PROBLEM_DESCRIPTION, NULL, mi->problemDescription, true);
LibFront::output_scalar(out, IS, 1, EMF_MATC_FILE_EMF, NULL, mi->matcInputFile_emf, true);
LibFront::output_scalar(out, IS, 1, EMF_MATC_FILE_SIF, NULL, mi->matcInputFile_sif, true);
//--If "save in model file"-flag is turned on (in Gui)
if ( mi->includePath_save )
LibFront::output_scalar(out, IS, 1, EMF_INCLUDE_PATH, NULL, mi->includePath, true);
if ( mi->resultsDirectory_save )
LibFront::output_scalar(out, IS, 1, EMF_RESULTS_DIRECTORY, NULL, mi->resultsDirectory, true);
if ( mi->temporaryFilesDirectory_save )
LibFront::output_scalar(out, IS, 1, EMF_LOG_DIRECTORY, NULL, mi->temporaryFilesDirectory, true);
LibFront::output_scalar(out, IS, 1, EMF_NOF_PROCESSORS, NULL, pi->nofProcessors);
LibFront::output_scalar(out, IS, 1, EMF_DIMENSION, NULL, mi->dimension);
LibFront::output_scalar(out, IS, 1, EMF_MINIMUM_EDGE_SIZE, NULL, mi->minEdgeSize);
if ( model->getFlagValue(GEOMETRY_EDITED_BODIES) ) {
LibFront::output_scalar(out, IS, 1, EMF_BODY_GEOMETRY_EDITED, NULL, 1);
}
if ( model->getFlagValue(GEOMETRY_EDITED_BOUNDARIES) ) {
LibFront::output_scalar(out, IS, 1, EMF_BOUNDARY_GEOMETRY_EDITED, NULL, 1);
}
//--Mesh names. These define targets for body grid paramters, boundary mesh-values etc.!
if ( mi->nofMeshes > 0 ) {
LibFront::output_vector(out, IS, 1, EMF_MESH_NAMES, "String", mi->nofMeshes, (const char**)mi->meshNames, true, true);
}
LibFront::output_scalar(out, IS, 1, EMF_CURRENT_MESH_INDEX, NULL, mi->currentMeshIndex);
//--Mesh control values and background mesh stuff ('mesh lists')
if ( cad_src && mi->nofMeshes > 0 ) {
LibFront::output_vector(out, IS, 1, EMF_MESH_H, NULL, mi->nofMeshes, mi->meshHs, false);
LibFront::output_vector(out, IS, 1, EMF_MESH_F, NULL, mi->nofMeshes, mi->meshFs, false);
}
//--Mesh background mesh stuff
if ( mi->nofBgMeshFiles > 0 ) {
LibFront::output_vector(out, IS, 1, EMF_MESH_BG_MESH_FILE_INDICES, NULL, mi->nofBgMeshFiles, mi->meshBgMeshFileIndices, false);
LibFront::output_vector(out, IS, 1, EMF_MESH_BG_MESH_FILES, "String", mi->nofBgMeshFiles, (const char**)mi->meshBgMeshFiles, true, true);
LibFront::output_vector(out, IS, 1, EMF_MESH_BG_MESH_ACTIVES, NULL, mi->nofBgMeshFiles, mi->meshBgMeshActives, false);
LibFront::output_vector(out, IS, 1, EMF_MESH_BG_MESH_CONTROLS, NULL, mi->nofBgMeshFiles, mi->meshBgMeshControls, false);
}
emf_outputSectionEnd(out);
return out;
}
// Output standard parameters.
ostream&
ModelOutputManager::emf_outputParameters(ostream& out, ecif_parameterType param_type)
{
const ModelInfo* mi = model->getModelInfo();
int index = 0;
while (true) {
Parameter* param = model->getParameter(index++, param_type);
if (param==NULL) break;
// Output
emf_outputSectionStart(out);
param->output_emf(out, ESF_INDENT_SIZE, 0, true, true, false);
emf_outputSectionEnd(out);
// Reset
param->resetValue();
}
return out;
}
#if 0
// Output all model vertex-points.
ostream&
ModelOutputManager::emf_outputVertices(ostream& out)
{
if (modelStatistics->nofVertices == 0)
return out;
short is = ESF_INDENT_SIZE;
BodyElement* v;
emf_outputSectionStart(out);
LibFront::output_string(out, is, 0, "Vertices");
// Variable: Vertex id
// Data: Coordinate point
int dim = 3; // Always 3D
int data_size = dim + 2; // Two ids; Bndr tag + BC-id
//LibFront::getFieldInfo(EFN_POINTS_AND_CONSTRAINT_IDS)->output_name(out, is, 1) << endl;
LibFront::indent(out, is, 1) << "Ids And Points" << endl;
LibFront::indent(out, is, 2) << "Variable Index" << endl;
LibFront::indent(out, is, 3) << "Size " << data_size << endl;
LibFront::indent(out, is, 4) << "Real" << endl;
v = getFirstVertex();
while ( v != NULL ) {
LibFront::indent(out, is, 5);
// Vertex tag
out << v->Tag() << " ";
// Boundary tag
out << v->getBoundaryTag() << " ";
// Boundary Condition Id
out << v->getBoundaryConditionId() << " ";
// Geometric point data
GcPoint* point = (GcPoint*)v->getGeometry();
out << *point;
out << endl;
v = getNextVertex();
}
// Data end marker
LibFront::indent(out, is, 4) << "End" << endl;
emf_outputSectionEnd(out);
return out;
}
#endif
// Data ouput after a section
ostream&
ModelOutputManager::emf_outputSectionEnd(ostream& out)
{
out << "End" << endl << endl;
return out;
}
// Data ouput before a section
ostream&
ModelOutputManager::emf_outputSectionStart(ostream& out)
{
//Currently nothing!
return out;
}
// Put statistics data into model output file.
ostream&
ModelOutputManager::emf_outputStatistics(ostream& out)
{
short is = ESF_INDENT_SIZE;
const ModelStatistics* ms = model->getModelStatistics();
emf_outputSectionStart(out);
// Header
LibFront::output_string(out, is, 0, "Statistics");
LibFront::output_scalar(out, is, 1, "Nof Bodies", NULL, ms->nofBodies);
LibFront::output_scalar(out, is, 1, "Nof Loops", NULL, ms->nofElementLoops);
// Calc nof body_elements to be output
int nof_elements = 0;
int index = 0;
while (true) {
BodyElement* be = model->getBoundary(index++);
if (be==NULL) break;
//-Check status
beStatus be_stat = be->getStatus();
if ( !(be_stat & (BE_DEVIDED | BE_SWAPPED)) )
nof_elements++;
}
LibFront::output_scalar(out, is, 1, "Nof Elements", NULL, nof_elements);
LibFront::output_scalar(out, is, 1, "Nof Outer Boundaries", NULL, ms->nofOuterBoundaries);
LibFront::output_scalar(out, is, 1, "Nof Inner Boundaries", NULL, ms->nofInnerBoundaries);
LibFront::output_scalar(out, is, 1, "Nof Vertices", NULL, ms->nofVertices);
// Calc max loop size
int max_loop_size = 0;
index = 0;
while (true) {
BodyElementLoop* bel = model->getBodyElementLoop(index++);
if (bel==NULL) break;
int size = bel->getNofElements();
if (size > max_loop_size)
max_loop_size = size;
}
LibFront::output_scalar(out, is, 1, "Max Loop Count", NULL, max_loop_size);
emf_outputSectionEnd(out);
return out;
}
// Put model timestamps data into ouput file.
ostream&
ModelOutputManager::emf_outputTimestamps(ostream& out)
{
char* QM = "\""; // quote-mark
short IS = ESF_INDENT_SIZE;
const ModelInfo* mi = model->getModelInfo();
emf_outputSectionStart(out);
// Timestamps when data was changed in the model file
// NOTE: Not necessarily the time when the target was updated!
// Header
LibFront::output_string(out, IS, 0, "Timestamps");
LibFront::output_scalar(out, IS, 1, "Front", NULL, mi->frontTs, true);
LibFront::output_scalar(out, IS, 1, "Database", NULL, mi->databaseTs, true);
LibFront::output_scalar(out, IS, 1, "Grid Parameter", NULL, mi->meshParameterTs, true);
LibFront::output_scalar(out, IS, 1, "Mesh", NULL, mi->meshTs, true);
LibFront::output_scalar(out, IS, 1, "Solver", NULL, mi->solverTs, true);
LibFront::output_scalar(out, IS, 1, "GebhardtFactors", NULL, mi->gebhardtFactorsTs, true);
LibFront::output_scalar(out, IS, 1, "ViewFactors", NULL, mi->viewfactorsTs, true);
emf_outputSectionEnd(out);
return out;
}
// Output model's vertex-table
//
ostream&
ModelOutputManager::emf_outputVertexTable(ostream& out)
{
VertexTable* vt = (VertexTable*)model->getVertexTable();
if ( vt == NULL || vt->dim1 * vt->dim2 == 0 ) return out;
emf_outputSectionStart(out);
LibFront::output_string(out, ESF_INDENT_SIZE, 0, EMF_VERTEX_TABLE, true);
LibFront::output_string(out, ESF_INDENT_SIZE, 1, EMF_POINTS, true);
LibFront::output_string(out, ESF_INDENT_SIZE, 2, EMF_SIZE, false);
out << " " << vt->dim1 << " " << vt->dim2 << endl;
LibFront::output_string(out, ESF_INDENT_SIZE, 3, EMF_REAL, true);
const char* def = getMatcString(vt->matcTable, EMF_POINTS);
//--Points defined by a Matc-function/variable
if ( model->keepMatcDefinitions() && def != NULL ) {
LibFront::output_matcDef(out, ESF_INDENT_SIZE, 3, NULL, NULL, def);
//--List of points
} else {
for ( int i = 0; i < vt->dim1; i++) {
BodyElement* v = model->getBodyElementById(vt->vertexIds[i]);
GcPoint* p = (GcPoint*)v->getGeometry();
LibFront::output_vector(out, ESF_INDENT_SIZE, 3, NULL, NULL, vt->dim2, (double*)p->getPoint(), false);
}
}
emf_outputSectionEnd(out);
return out;
}
void
ModelOutputManager::initClass(Model* mdl)
{
ModelOutputManager::model = mdl;
}
//******************************************
//*** MESH INPUT format output functions ***
//******************************************
// Form Mesh Input File (mif-file)
//
ostream&
ModelOutputManager::mif_output(ostream& out)
{
const ModelInfo* mi = model->getModelInfo();
const ModelStatistics* ms = model->getModelStatistics();
int i;
//--Set tags for possible linearizing boundary points
// NOTE: also mark all points inactive for mesh
// Point to be included are marked active below via bodies!
//--Also set possible boundary point mesh density values
//
model->setBoundaryPointData();
int mesh_index = mi->currentMeshIndex;
//--Find active mesh objects
bool* active_objects = new bool[1 + ms->nofModelObjects];
// Init flags
for (i = 0; i <= ms->nofModelObjects; i++) {
active_objects[i] = false;
}
int index = 0;
// Mark active object via bodies!
while (true) {
Body* body = model->getBody(index++);
if ( body == NULL ) break;
if ( body->isExcludedFromMesh(mesh_index) ) continue;
body->markActiveMeshObjects(mesh_index, active_objects);
}
//--Count active objects by type
int nof_vertices = 0;
int nof_edges = 0;
int nof_bodies = 0;
int nof_inner_boundaries = 0;
int nof_outer_boundaries = 0;
for (i = 0; i <= ms->nofModelObjects; i++) {
if ( !active_objects[i] )
continue;
Body* body = NULL;
BodyElement* be = NULL;
ModelObject* obj = model->getModelObjectById(i);
int bd1_id, bd2_id;
if ( obj == NULL )
continue;
switch ( obj->getObjectType() ) {
case OT_VERTEX:
nof_vertices++;
break;
case OT_EDGE:
be = model->getBodyElementById(i);
nof_edges += be->getNofMifGeometries();
break;
case OT_BODY:
body = model->getBodyById(i);
if ( !body->isVirtual() && !body->isExcludedFromMesh(mesh_index) ) {
nof_bodies++;
}
break;
}
}
int nof_bpoints = mif_getNofBoundaryPoints();
//--Output file
// Comments
// ========
char tb[128];
getCurrentTs(tb, 128);
out << "!ElmerMesh input file from ElmerFront" << endl;
out << "!Saved = " << tb << endl;
out << "!Case = " << mi->modelName << " " << mi->problemName << endl;
out << "!Model dir = " << mi->modelDirectory_absolute << endl;
out << "!Include path = " << mi->includePath_absolute << endl;
out << "!Results dir = " << mi->resultsDirectory_absolute << endl;
double mesh_h = mi->meshHs[mesh_index];
double mesh_f = mi->meshFs[mesh_index];
// Header
// ======
out << "Geometry2D:" << endl;
if ( mesh_h > 0 )
indent(out, 2) << "H: " << mesh_h << endl;
if ( mesh_f > 0 )
indent(out, 2) << "MeshScalingFactor: " << mesh_f << endl;
indent(out, 2) << "Nodes: " << nof_vertices + nof_bpoints << endl;
indent(out, 2) << "Edges: " << nof_edges << endl;
indent(out, 2) << "Bodies: " << nof_bodies << endl;
// Vertices + BoundaryPoints
// =========================
//out << "!Vertices: id, x, y, meshH" << endl;
//-Vertices
index = 0;
while (true) {
BodyElement* v = model->getVertex(index++);
if (v==NULL) break;
if ( active_objects[v->Id()] ) {
v->output_mif(out);
}
}
//-BoundaryPoints
mif_outputBoundaryPoints(out);
// Edges
// =====
//out << "!Edges: id, meshH, meshN, vertex count, vertex ids" << endl;
index = 0;
int next_mtag = 1;
while (true) {
BodyElement* e = model->getEdge(index++);
if (e==NULL) break;
if ( active_objects[e->Id()] ) {
e->setMifTag(next_mtag);
e->output_mif(out);
}
}
// Bodies
// ======
index = 0;
while (true) {
Body* b = model->getBody(index++);
if ( b == NULL ) break;
if ( b->isExcludedFromMesh(mesh_index) ) continue;
if ( !b->isVirtual() && active_objects[b->Id()] ) {
b->output_mif(out);
}
}
// End of file
// ===========
out << "!End" << endl;
return out;
}
// Get nof active mesh boundary points
//
int
ModelOutputManager::mif_getNofBoundaryPoints()
{
int i, index;
int active_count = 0;
int bp_count;
BoundaryPoint** bp_points;
//--First mark all unchecked (to handle copied points)
//
index = 0;
while (true) {
BodyElement* be = model->getBoundary(index++);
if (be==NULL) break;
bp_count = 0;
bp_points = NULL;
be->getBoundaryPoints(bp_count, bp_points);
for (i = 0; i < bp_count; i++) {
BoundaryPoint* bp = bp_points[i];
// If this is an original vertex, we do not mark it
if ( bp->isVertex() ) {
continue;
}
bp->checked = false;
}
delete[] bp_points;
bp_points = NULL;
}
//--Next count nof actives in meshing
//
index = 0;
while (true) {
BodyElement* be = model->getBoundary(index++);
if (be==NULL) break;
bp_count = 0;
bp_points = NULL;
be->getBoundaryPoints(bp_count, bp_points);
for (i = 0; i < bp_count; i++) {
BoundaryPoint* bp = bp_points[i];
if ( !bp->activeInMeshing ||
bp->isVertex() ||
bp->checked
) {
continue;
}
active_count++;
// Mark as handled
bp->checked = true;
}
delete[] bp_points;
bp_points = NULL;
}
return active_count;
}
// Output linearizing boundary points (as nodes) to
// mesh input file (mif-file)
// NOTE: original vertices are not output here, they
// are output separately (before these points)
//
ostream&
ModelOutputManager::mif_outputBoundaryPoints(ostream& out)
{
int i, index;
Point3 p;
bool only_active = true;
// For mesh density info
int mesh_index = model->getCurrentMeshIndex();
int bp_count;
BoundaryPoint** bp_points;
//--First mark all unchecked to handle copied (ex. loop end) points
//
index = 0;
while (true) {
BodyElement* be = model->getBoundary(index++);
if (be==NULL) break;
bp_count = 0;
bp_points = NULL;
be->getBoundaryPoints(bp_count, bp_points);
for (i = 0; i < bp_count; i++) {
BoundaryPoint* bp = bp_points[i];
// If this is an original vertex, we do not mark
if ( bp->isVertex() ) {
continue;
}
bp->checked = false;
}
delete[] bp_points;
bp_points = NULL;
}
//---Next output all active
//
index = 0;
while (true) {
BodyElement* be = model->getBoundary(index++, only_active);
if (be==NULL) break;
bp_count = 0;
bp_points = NULL;
be->getBoundaryPoints(bp_count, bp_points);
for (i= 0; i < bp_count; i++) {
BoundaryPoint* bp = bp_points[i];
if ( !bp->activeInMeshing ||
bp->isVertex() ||
bp->checked
) {
continue;
}
bp->point->getPoint(p);
bp->checked = true;
out << "NodeId: " << bp->tag << " -1";
if ( bp->meshDensityType != ' ' ) {
out << " " << bp->meshDensityType << ": " << bp->meshDensityValue;
}
// NOTE: Only two coordinates are output for Mesh2D
//
out << " " << p[0] << " " << p[1];
out << endl;
}
delete[] bp_points;
bp_points = NULL;
}
return out;
}
//*******************************
//*** MATC definitions output ***
//*******************************
ostream&
ModelOutputManager::outputMatcDefinitions(ostream& out, int nof_defs, char** defs, bool dsign)
{
int i;
// First write variables
for (i = 0; i < nof_defs; i++) {
if ( 0 == strncmp("function", LibFront::trimLeft(defs[i]), 8) ) continue;
if ( dsign ) out << '$';
out << defs[i] << endl;
}
// Then write functions
for (i = 0; i < nof_defs; i++) {
if ( 0 != strncmp("function", LibFront::trimLeft(defs[i]), 8) ) continue;
if ( dsign ) out << '$';
out << defs[i] << endl;
}
return out;
}
//********************************************
//*** SOLVER INPUT format output functions ***
//********************************************
// Form Solver Input File (sif-file)
//
ostream&
ModelOutputManager::sif_output(ostream& out)
{
UserInterface* gui = theControlCenter->getGui();
Parameter* p;
ParameterField* pf;
int index;
const ModelInfo* mi = model->getModelInfo();
const ModelStatistics* ms = model->getModelStatistics();
//-Calculate nof active solvers
int nof_solvers = 0;
index = 0;
while (true) {
p = model->getParameter(index++, ECIF_SOLVER);
if (p==NULL) break;
pf = p->getFieldBySifName("Active");
if ( p->IsActive() ) {
nof_solvers++;
}
}
//-Calculate nof active calculators
int nof_calculators = 0;
index = 0;
while (p != NULL) {
p = model->getParameter(index++, ECIF_CALCULATOR);
if (p==NULL) break;
pf = p->getFieldBySifName("Active");
if ( p->IsActive() ) {
nof_calculators++;
}
}
// Info as comments
out << "!ElmerSolver input file from ElmerFront" << endl;
out << "!Saved = " << mi->modified << endl;
out << "!Case = " << mi->modelName << " " << mi->problemName << endl;
out << "!Model dir = " << mi->modelDirectory_absolute << endl;
out << "!Include path = " << mi->includePath_absolute << endl;
out << "!Results dir = " << mi->resultsDirectory_absolute << endl;
out << endl;
// NOTE: output "max counters" for Solver, ie highest id in use!
// NOT: These are not used any more by solver (Elmer3.0-->), so printed as comments
//
out << "!Bodies " << ms->nofBodies << endl;
out << "!Equations " << ms->nofEquations << endl;
out << "!Solvers " << nof_solvers + nof_calculators << endl;
out << "!Materials " << ms->nofMaterials << endl;
out << "!Body Forces " << ms->nofBodyForces << endl;
out << "!Initial Conditions " << ms->nofInitialConditions << endl;
out << "!Boundary Conditions " << ms->nofBoundaryConditions << endl;
out << "!Boundaries " << ms->nofInnerBoundaries + ms->nofOuterBoundaries << endl;
out << endl;
// Sif file echo on/off
p = model->getParameterById(ECIF_SOLVER_CONTROL, 1);
if ( p != NULL ) {
bool echo_on = false;
p->getFieldValueBySifName(SIF_ECHO_ON, echo_on);
// If not on, make it a comment!
if ( !echo_on ) {
out << "!";
}
out << "echo on" << endl;
out << endl;
}
sif_outputHeader(out);
// Possible matc input file
//
if ( mi->matcInputFile_sif != NULL &&
mi->matcInputFile_sif[0] != '\0'
) {
out << "!Matc input file" << endl;
out << "$source(\"" << mi->matcInputFile_sif << "\")" << endl;
out << endl;
}
// Possible matc definitions
//
int nof_defs = 0;
char** defs = NULL;
gui->getMatcSifDefinitions(nof_defs, defs);
// Write matc definitions into Sif-file (true <--> with $-sign)
if ( nof_defs > 0 ) {
out << "!Matc definitions" << endl;
outputMatcDefinitions(out, nof_defs, defs, true);
out << endl;
}
// Possible model parameter
//
p = model->getParameterById(ECIF_MODEL_PARAMETER, 1);
if ( p != NULL ) {
SifOutputControl soc;
soc.outputId = false;
soc.outputName = false;
soc.outputType = false;
soc.outputAll = false;
soc.sectionName = SIF_HEADER;
p->output_sif(out, ESF_INDENT_SIZE, 0, soc);
out << endl;
}
sif_outputSimulation(out);
sif_outputConstants(out);
sif_outputBodies(out);
sif_outputEquations(out);
sif_outputSolvers(out);
sif_outputMaterials(out);
sif_outputBodyForces(out);
sif_outputInitialConditions(out);
sif_outputBoundaryConditions(out);
sif_outputBoundaries(out);
sif_outputEof(out);
return out;
}
//Data after block header
ostream&
ModelOutputManager::sif_outputAfterHeader(ostream& out)
{
out << endl << "!" << endl;
return out;
}
// Output bodies.
ostream&
ModelOutputManager::sif_outputBodies(ostream& out)
{
int index = 0;
while (true) {
Body* body = model->getBody(index++);
if (body==NULL) break;
sif_outputSectionStart(out);
body->output_sif(out, ESF_INDENT_SIZE, 0);
sif_outputSectionEnd(out);
}
return out;
}
//----- ESF Output functions (Elemer Solver input format)
// Output body-forces data.
ostream&
ModelOutputManager::sif_outputBodyForces(ostream& out)
{
SifOutputControl soc;
soc.sectionName = SIF_BODY_FORCE;
int index = 0;
while (true) {
Parameter* force = model->getParameter(index++, ECIF_BODY_FORCE);
if (force==NULL) break;
if ( force->getApplyCount() > 0 ) {
sif_outputSectionStart(out);
force->output_sif(out, ESF_INDENT_SIZE, 0, soc);
sif_outputSectionEnd(out);
}
}
return out;
}
// output All Boundaries.
ostream&
ModelOutputManager::sif_outputBoundaries(ostream& out)
{
//BodyElement* be = model->getFirstBodyElement();
int index = 0;
while (true) {
BodyElement* be = model->getBoundary(index++);
if (be==NULL) break;
//-Check status
beStatus be_stat = be->getStatus();
if ( be_stat & (BE_DEVIDED | BE_SWAPPED) ) {
continue;
}
// NOTE: Only those boundaries are output which have
// a boundary parameter defined (this is the only way
// to output boundary parameters!)
//
if ( NO_INDEX == be->getBoundaryParameterId() ) {
continue;
}
sif_outputSectionStart(out);
be->output_sif(out, ESF_INDENT_SIZE, 0);
sif_outputSectionEnd(out);
}
return out;
}
#if 0
// output All Boundaries.
ostream&
ModelOutputManager::sif_outputBoundaries(ostream& out)
{
int count, i;
int boundary_condition_id;
// Outer boundaries.
OuterBoundary* ob;
count = modelData->outerBoundaries->size();
for (i = 0; i < count; i++) {
ob = (*modelData->outerBoundaries)[i];
if (ob == 0)
continue;
boundary_condition_id = ob->outerBoundary->getConditionId();
sif_outputSectionStart(out);
out << "Boundary " << ob->outerBoundary->Tag()
<< endl;
out << LibFront::indent(ESF_INDENT_SIZE, 1) << "Body 1"
<< endl;
out << LibFront::indent(ESF_INDENT_SIZE, 2) << "Integer "
<< ob->parentBody->Tag()
<< endl;
// NOTE: bc-cond-id -1 is also output!
out << endl;
out << LibFront::indent(ESF_INDENT_SIZE, 1) << "Boundary Condition"
<< endl;
out << LibFront::indent(ESF_INDENT_SIZE, 2) << "Integer "
<< boundary_condition_id
<< endl;
sif_outputSectionEnd(out);
}
// Inner boundaries
InnerBoundary* ib;
count = modelData->innerBoundaries->size();
for (i = 0; i < count; i++) {
ib = (*modelData->innerBoundaries)[i];
if (ib == 0)
continue;
boundary_condition_id = ib->innerBoundary->getConditionId();
sif_outputSectionStart(out);
out << "Boundary " << ib->innerBoundary->Tag()
<< endl;
out << LibFront::indent(ESF_INDENT_SIZE, 1) << "Body 1"
<< endl;
out << LibFront::indent(ESF_INDENT_SIZE, 2) << "Integer "
<< ib->parentPair[0]->body->Tag()
<< endl;
out << endl;
out << LibFront::indent(ESF_INDENT_SIZE, 1) << "Body 2"
<< endl;
out << LibFront::indent(ESF_INDENT_SIZE, 2) << "Integer "
<< ib->parentPair[1]->body->Tag()
<< endl;
if (boundary_condition_id != NO_INDEX) {
out << endl;
out << LibFront::indent(ESF_INDENT_SIZE, 1) << "Boundary Condition"
<< endl;
out << LibFront::indent(ESF_INDENT_SIZE, 2) << "Integer "
<< boundary_condition_id
<< endl;
}
sif_outputSectionEnd(out);
}
return out;
}
#endif
// Output all boundary conditions (outer, inner).
ostream&
ModelOutputManager::sif_outputBoundaryConditions(ostream& out)
{
SifOutputControl soc;
soc.sectionName = SIF_BOUNDARY_CONDITION;
int index = 0;
while (true) {
Parameter* cond = model->getParameter(index++, ECIF_BOUNDARY_CONDITION);
if (cond==NULL) break;
if (cond->getApplyCount() > 0 ) {
cond->updateTargetTags();
sif_outputSectionStart(out);
cond->output_sif(out, ESF_INDENT_SIZE, 0, soc);
sif_outputSectionEnd(out);
}
}
return out;
}
// Output physial constants.
ostream&
ModelOutputManager::sif_outputConstants(ostream& out)
{
sif_outputSectionCode(out, SIF_CONSTANT);
SifOutputControl soc;
soc.outputType = false;
soc.outputName = false;
soc.sectionName = SIF_CONSTANT;
int index =0;
while (true) {
Parameter* constant = model->getParameter(index++, ECIF_CONSTANT);
if (constant==NULL) break;
sif_outputSectionStart(out);
constant->output_sif(out, ESF_INDENT_SIZE, 0, soc);
sif_outputSectionEnd(out);
}
return out;
}
// Output all initial conditions.
ostream&
ModelOutputManager::sif_outputInitialConditions(ostream& out)
{
SifOutputControl soc;
soc.sectionName = SIF_INITIAL_CONDITION;
int index = 0;
while (true) {
Parameter* cond = model->getParameter(index++, ECIF_INITIAL_CONDITION);
if (cond==NULL) break;
if (cond->getApplyCount() > 0 ) {
sif_outputSectionStart(out);
cond->output_sif(out, ESF_INDENT_SIZE, 0, soc);
sif_outputSectionEnd(out);
}
}
return out;
}
// Mark End Of File.
ostream&
ModelOutputManager::sif_outputEof(ostream& out)
{
out << "!End Of File" << endl;
return out;
}
// Output equations data.
ostream&
ModelOutputManager::sif_outputEquations(ostream& out)
{
SifOutputControl soc;
soc.sectionName = SIF_EQUATION;
int index = 0;
while (true) {
Parameter* equation = model->getParameter(index++, ECIF_EQUATION);
if (equation==NULL) break;
if ( equation->getApplyCount() > 0 ) {
sif_outputSectionStart(out);
equation->output_sif(out, ESF_INDENT_SIZE, 0, soc);
sif_outputSectionEnd(out);
}
}
return out;
}
// Put model header data into ouput file.
ostream&
ModelOutputManager::sif_outputHeader(ostream& out)
{
Parameter* p;
ParameterField* pf;
char key_buffer[1024];
char value_buffer[1024];
char QM = '\"';
const ModelInfo* mi = model->getModelInfo();
sif_outputSectionStart(out);
sif_outputSectionCode(out, SIF_HEADER);
// Keyword checking (ignore, warn, abort)
p = model->getParameterById(ECIF_SOLVER_CONTROL, 1);
if ( p != NULL ) {
LibFront::toUpper(SIF_CHECK_KEYWORDS, key_buffer);
p->getFieldValueBySifName(SIF_CHECK_KEYWORDS, 1023, value_buffer);
indent(out, ESF_INDENT_SIZE);
out << key_buffer << " " << value_buffer;
out << endl;
}
// Meshdir
LibFront::indent(out, ESF_INDENT_SIZE, 1)
<< "Mesh DB "
<< QM << model->getMeshDirValue() << QM // "MESHDIR"
<< " ";
if ( mi->nofActiveMeshes > 0 &&
mi->meshNames != NULL
) {
out << QM << mi->meshNames[mi->activeMeshIndices[0]] << QM;
} else {
out << QM << QM;
}
out << endl;
LibFront::indent(out, ESF_INDENT_SIZE, 1) << "Include Path " << QM << mi->includePath << QM << endl;
LibFront::indent(out, ESF_INDENT_SIZE, 1) << "Results Directory " << QM << mi->resultsDirectory << QM << endl;
sif_outputSectionEnd(out);
return out;
}
// Output materials data.
ostream&
ModelOutputManager::sif_outputMaterials(ostream& out)
{
SifOutputControl soc;
soc.sectionName = SIF_MATERIAL;
int index = 0;
while (true) {
Parameter* mater = model->getParameter(index++, ECIF_MATERIAL);
if (mater==NULL) break;
if ( mater->getApplyCount() ) {
sif_outputSectionStart(out);
mater->output_sif(out, ESF_INDENT_SIZE, 0, soc);
sif_outputSectionEnd(out);
}
}
return out;
}
// Output code for a new section (SEC_HEADER etc.)
ostream&
ModelOutputManager::sif_outputSectionCode(ostream& out, const char* section_cd)
{
out << section_cd << endl;
return out;
}
// Data ouput after a section
ostream&
ModelOutputManager::sif_outputSectionEnd(ostream& out)
{
out << "End" << endl << endl;
return out;
}
// Data ouput before a section
ostream&
ModelOutputManager::sif_outputSectionStart(ostream& out)
{
//Currently nothing!
return out;
}
// Put problem data into ouput file.
ostream&
ModelOutputManager::sif_outputSimulation(ostream& out)
{
Parameter* p;
const ModelInfo* mi = model->getModelInfo();
sif_outputSectionStart(out);
sif_outputSectionCode(out, SIF_SIMULATION);
SifOutputControl soc;
soc.outputType = false;
soc.outputName = false;
soc.sectionName = SIF_SIMULATION;
// Solver control
p = model->getParameterById(ECIF_SOLVER_CONTROL, 1);
if ( p != NULL ) {
SifOutputControl soc;
soc.outputId = false;
soc.outputName = false;
soc.outputType = false;
soc.outputAll = false;
soc.sectionName = SIF_SIMULATION;
p->output_sif(out, ESF_INDENT_SIZE, 0, soc);
}
// Simulation parameter
p = model->getParameterById(ECIF_SIMULATION_PARAMETER, 1);
if ( p != NULL ) {
SifOutputControl soc;
soc.outputId = false;
soc.outputName = false;
soc.outputType = false;
soc.outputAll = false;
soc.sectionName = SIF_SIMULATION;
p->output_sif(out, ESF_INDENT_SIZE, 0, soc);
out << endl;
}
// Coordinate stuff
Parameter* coordinate = model->getParameter(0, ECIF_COORDINATE);
if (coordinate != NULL ) {
sif_outputSectionStart(out);
coordinate->output_sif(out, ESF_INDENT_SIZE, 0, soc);
out << endl;
}
// Timestep stuff (first active!)
int index = 0;
while (true) {
Parameter* timestep = model->getParameter(index++, ECIF_TIMESTEP);
if (timestep==NULL) break;
ParameterField* pf = timestep->getFieldBySifName("Active");
// If active, output and stop looping
if ( timestep->IsActive()) {
sif_outputSectionStart(out);
timestep->output_sif(out, ESF_INDENT_SIZE, 1, soc);
break;
}
}
out << endl;
// Result and inputfile stuff
Parameter* datafile = model->getParameter(0, ECIF_DATAFILE);
if ( datafile != NULL ) {
sif_outputSectionStart(out);
datafile->output_sif(out, ESF_INDENT_SIZE, 0, soc);
}
// Mesh Input File name
char* mif_file_name = NULL;
model->getMeshInputFileName(mif_file_name, 0);
if ( mif_file_name != NULL ) {
LibFront::output_scalar(out, ESF_INDENT_SIZE, 1, "Mesh Input File", "File", mif_file_name, true);
out << endl;
delete[] mif_file_name;
}
sif_outputSectionEnd(out);
return out;
}
// Output solvers data.
ostream&
ModelOutputManager::sif_outputSolvers(ostream& out)
{
int i;
int index;
Parameter* solver;
Parameter* calculator;
const ModelStatistics* ms = model->getModelStatistics();
int count = ms->nofSolvers + ms->nofCalculators;
Parameter** solvers = new Parameter*[count];
for (i = 0; i < count; i++) {
solvers[i] = NULL;
}
// Set active solvers and calculators according to the SOLVING_ORDER field
// into the solvers table for output
int counter = 0;
//-Solvers
index = 0;
while (true) {
solver = model->getParameter(index++, ECIF_SOLVER);
if (solver==NULL) break;
ParameterField* order_fld = solver->getFieldBySifName("Solving Order");
if ( solver->IsActive() ) {
counter++;
int index;
if (order_fld != NULL)
index = atol(order_fld->getDataStrings()[0]);
else
index = counter;
solvers[index - 1] = solver;
}
}
//-Calculators
index = 0;
while (true) {
calculator = model->getParameter(index++, ECIF_CALCULATOR);
if (calculator==NULL) break;
if ( calculator->IsActive() ) {
counter++;
int index;
ParameterField* order_fld = calculator->getFieldBySifName("Solving Order");
if (order_fld != NULL)
index = atol(order_fld->getDataStrings()[0]);
else
index = counter;
solvers[index - 1] = calculator;
}
}
SifOutputControl soc;
soc.outputName = false;
soc.sectionName = SIF_SOLVER;
Parameter* sc = model->getParameter(0, ECIF_SOLVER_CONTROL);
ParameterField* pf = NULL;
if ( sc != NULL ) {
pf = sc->getFieldByGuiName("RELOAD_INPUT_FILE", true);
}
// If Reload-file field active, output reload info as a Solver
//
if ( pf != NULL ) {
soc.reloadSolverIsOutput = true;
LibFront::output_string(out, ESF_INDENT_SIZE, 0, "Solver 1", true);
pf = sc->getFieldByGuiName("EXEC_SOLVER", true);
if ( pf != NULL ) {
pf->output_sif(out, ESF_INDENT_SIZE, 1, "Solver", true);
}
// Equation for the 'solver'
if ( model->getSolverKeywordTypeGiven("Solver", "Equation" ) )
LibFront::output_string(out, ESF_INDENT_SIZE, 1, "Equation = \"ReloadInput\"", true);
else
LibFront::output_string(out, ESF_INDENT_SIZE, 1, "Equation = String \"ReloadInput\"", true);
// Procedure for the 'solver'
if ( model->getSolverKeywordTypeGiven("Solver", "Procedure" ) )
LibFront::output_string(out, ESF_INDENT_SIZE, 1, "Procedure = \"ReloadInput\" \"ReloadInput\"", true);
else
LibFront::output_string(out, ESF_INDENT_SIZE, 1, "Procedure = File \"ReloadInput\" \"ReloadInput\"", true);
LibFront::output_string(out, ESF_INDENT_SIZE, 0, "End", true);
out << endl;
}
// Output active solvers in solving order
for (i = 0; i < count; i++) {
if ( solvers[i] == NULL )
continue;
sif_outputSectionStart(out);
solvers[i]->output_sif(out, ESF_INDENT_SIZE, 0, soc);
sif_outputSectionEnd(out);
}
delete[] solvers;
return out;
}
// Output those equation fields in whose target is actually a solver section
//
ostream&
ModelOutputManager::sif_outputSolverTargetFields(ostream& out, short indent_size, short indent_level, const char* source_eq_name)
{
NameSet targetFieldNames;
int nof_values;
bool* values = NULL;
int index = 0;
while (true) {
delete[] values;
values = NULL;
Parameter* equation = model->getParameter(index++, ECIF_EQUATION);
if (equation==NULL) break;
if ( equation->getApplyCount() == 0 ) continue;
if ( !equation->getFieldValueBySifName(source_eq_name, nof_values, values) ) {
continue;
}
if ( !values[0] ) continue;
equation->outputSolverTargetFields_sif(out, indent_size, indent_level, source_eq_name, targetFieldNames);
}
return out;
}
//***********************************
//*** MESH FILES output functions ***
//***********************************
// Elmer format file
// =================
void
ModelOutputManager::write_Elmer_mesh(char* mesh_dir)
{
const ModelInfo* mi = model->getModelInfo();
const MeshInfo* mei = model->getMeshInfo();
const ModelStatistics* ms = model->getModelStatistics();
UserInterface* gui = theControlCenter->getGui();
int info;
int nodeC = mei->nofNodes;
int elemC_net = mei->nofBulkElements - mei->nofSplittedElements; // To header!
int belemC = mei->nofBoundaryElements;
int elemC = mei->nofBulkElements; // For looping bulk table!
int counter;
int i, index;
// These are 101-type boundary elements!
mei->nofUsedElementTypes[0] = ms->nofVertices;
// Count nof element types used
int typeC = 0;
for (i = 0; i < MAX_NOF_ELEM_CODES; i++) {
if ( mei->nofUsedElementTypes[i] == 0 )
continue;
typeC++;
}
int* types = new int[typeC];
int* typesC = new int[typeC];
// Count nof elements by used types
counter = 0;
for (i = 0; i < MAX_NOF_ELEM_CODES; i++) {
if ( mei->nofUsedElementTypes[i] == 0 )
continue;
types[counter] = MeshElementDesc[i][DESC_ELEM_TYPE];
typesC[counter] = mei->nofUsedElementTypes[i];
// Update boundary element counter
if ( mi->dimension == ECIF_3D ) {
// Add possible edges and vertices
if ( types[counter] < 300 ) {
belemC += typesC[counter];
}
} else {
// Add possible vertices
if ( types[counter] < 200 ) {
belemC += typesC[counter];
}
}
counter++;
}
// Init DB
eio_init(info);
// EIO (create mesh files)
// -----------------------
eio_create_mesh(mesh_dir, info);
// Ouch, no success with the mesh!
if (info == -1) {
gui->errMsg(0, "Cannot create mesh in database using directory and name:", mesh_dir);
return;
}
// EIO (mesh.header)
// -----------------
eio_set_mesh_description(nodeC, elemC_net, belemC, typeC, types, typesC, info);
// Nodes
// =====
Point3* nodes = model->getMeshNodeData();
double coordinates[3];
for (i = 0; i < nodeC; i++) {
for (short j = 0; j < 3; j++) {
coordinates[j] = nodes[i][j];
}
// Elmer DB ids
int node_id = i + 1;
int cond_id = 0;
// EIO (mesh.nodes)
// ----------------
eio_set_mesh_node(node_id, cond_id, coordinates, info);
}
int node_ids[32];
// Elements
// ========
// Bulk elements
// -------------
MeshElementTable* bt = model->getMeshBulkElements();
counter = 0;
for (i = 0; i < elemC; i++) {
if (bt->splitted[i]) {
continue;
}
int elem_code = bt->getElementCode(i);
int elem_type = MeshElementDesc[elem_code][DESC_ELEM_TYPE];
int body_id = model->getModelObjectTagById(bt->parentIds[i][0]);
int nof_nodes = MeshElementDesc[elem_code][DESC_NOF_NODES];
const int* elem_node_ids = bt->getNodeIds(i);
// Elmer DB ids
int elem_ext_id = 1 + counter++;
for (short j = 0; j < nof_nodes; j++) {
node_ids[j] = 1 + elem_node_ids[j];
}
// EIO (mesh.elements)
// -------------------
eio_set_mesh_element_conns(elem_ext_id, body_id, elem_type, node_ids, info);
}
// Boundary elements
// -----------------
counter = 0;
BodyElement* be;
//-Boundary elements proper (faces or edges)
MeshElementTable* bet = model->getMeshBoundaryElements();
index = 0;
while (true) {
be = model->getBoundary(index++);
if (be==NULL) break;
int nof_elements = be->getNofMeshElements();
for (i = 0; i < nof_elements; i++) {
int index = be->getMeshElementId(i);
// parent element ids, Elmer DB numbering
int p_elem1_int_id = bet->parentIds[index][0]; // Int parent elem1
int p_elem2_int_id = bet->parentIds[index][1]; // Int parent elem2
if ( model->modelHasBulkRenumbering() ) {
p_elem1_int_id = model->getRenumberedMeshBulkElementId(p_elem1_int_id);
p_elem2_int_id = model->getRenumberedMeshBulkElementId(p_elem2_int_id);
}
int p_elem1_ext_id = 1 + p_elem1_int_id; // Ext parent elem1 id
int p_elem2_ext_id = 1 + p_elem2_int_id; // Ext parent elem2 id
// Element code for the boundary element
meshElementCode elem_code = bet->getElementCode(index);
int elem_type = MeshElementDesc[elem_code][DESC_ELEM_TYPE];
int nof_nodes = MeshElementDesc[elem_code][DESC_NOF_NODES];
// Boundary tag for the element
int bndr_tag = be->getBoundaryTag(); // Internal boundary tag
const int* elem_node_ids = bet->getNodeIds(index);
// Elmer DB ids (DB index origo is 1)
for (int k = 0; k < nof_nodes; k++) {
node_ids[k] = 1 + (int)elem_node_ids[k];
}
//int bndr_elem_id = index + 1;
int elem_id = ++counter;
info = 0;
// EIO (mesh.boundary)
// -------------------
eio_set_mesh_bndry_element(elem_id, bndr_tag,
p_elem1_ext_id, p_elem2_ext_id,
elem_type, node_ids,
info);
} // for all mesh elements in the boundary
} // All boundaries
//-Boundary elements edges (3D models)
//
MeshElementTable* eet = model->getMeshBoundaryElementEdges();
index = 0;
while ( mi->dimension == ECIF_3D) {
be = model->getEdge(index++);
if (be==NULL) break;
int nof_elements = be->getNofMeshElements();
for (i = 0; i < nof_elements; i++) {
int index = be->getMeshElementId(i);
// No parents for edges in 3D!
int p_elem1_ext_id = -1;
int p_elem2_ext_id = -1;
// Element code for the boundary element
meshElementCode elem_code = eet->getElementCode(index);
int elem_type = MeshElementDesc[elem_code][DESC_ELEM_TYPE];
int nof_nodes = MeshElementDesc[elem_code][DESC_NOF_NODES];
// Boundary tag for the element
int bndr_tag = be->getBoundaryTag(); // Internal boundary tag
const int* elem_node_ids = eet->getNodeIds(index);
// Elmer DB ids (DB index origo is 1)
for (int k = 0; k < nof_nodes; k++) {
node_ids[k] = 1 + elem_node_ids[k];
}
//int bndr_elem_id = index + 1;
int elem_id = ++counter;
info = 0;
// EIO (mesh.boundary)
// -------------------
eio_set_mesh_bndry_element(elem_id, bndr_tag,
p_elem1_ext_id, p_elem2_ext_id,
elem_type, node_ids,
info);
} // for all meshelements in the edge
} // All 3D edges
// Vertices
// =================
index = 0;
while (true) {
be = model->getVertex(index++);
if (be==NULL) break;
int elem_id = ++counter;
int bndr_tag = be->getBoundaryTag();
int p_elem1_ext_id = -1;
int p_elem2_ext_id = -1;
int elem_type = 101;
node_ids[0] = 1 + be->getMeshElementId(0);
// EIO (mesh.boundary)
// -------------------
eio_set_mesh_bndry_element(elem_id, bndr_tag,
p_elem1_ext_id, p_elem2_ext_id,
elem_type, node_ids,
info);
}
// Close EIO
eio_close_mesh(info);
eio_close(info);
delete[] types;
delete[] typesC;
gui->showMsg("Mesh saved in Elmer format!");
}
// Elmer Post format file
// ======================
void
ModelOutputManager::write_ElmerPost_mesh(ostream& outfile)
{
int i, index;
const ModelInfo* mi = model->getModelInfo();
const MeshInfo* mei = model->getMeshInfo();
const ModelStatistics* ms = model->getModelStatistics();
UserInterface* gui = theControlCenter->getGui();
int nofNodes = mei->nofNodes;
int nofElements = mei->nofBulkElements - mei->nofSplittedElements;
int nofBoundaryElements = mei->nofBoundaryElements;
// Header
// =====
outfile << nofNodes << " " << nofElements + nofBoundaryElements;
outfile << endl;
// Node coordinates
// ================
Point3* nodes = model->getMeshNodeData();
for (i = 0; i < nofNodes; i++) {
Point3& p = nodes[i];
outfile << " " << p[0] << " " << p[1] << " " << p[2] << endl;
}
// Bulk elements from bodies
// =========================
MeshElementTable* bt = model->getMeshBulkElements();
outfile << "#group all-bodies" << endl;
while (true) {
Body* body = model->getBody(index++);
if (body==NULL) break;
int nof_elements = body->getNofMeshElements();
if (nof_elements == 0) {
continue;
}
//--Body name as the group name
outfile << "#group " << body->getName() << endl;
// All body mesh elements
for (i = 0; i < nof_elements; i++) {
int eid = body->getMeshElementId(i);
const int* nodes = bt->getNodeIds(eid);
meshElementCode elem_code = bt->getElementCode(eid);
int body_id = bt->parentIds[eid][0];
outfile << "all ";
// Elm type for the bulk element
outfile << MeshElementDesc[elem_code][DESC_ELEM_TYPE];
// Internal node ids
outfile << " ";
int nof_nodes = MeshElementDesc[elem_code][DESC_NOF_NODES];
for (int j = 0; j < nof_nodes; j++) {
outfile << " " << nodes[j];
}
outfile << endl;
}
//--End for this body group
outfile << "#endgroup " << body->getName() << endl;
}
//---End for all bodies
outfile << "#endgroup all-bodies" << endl;
// Boundary elements from boundaries
// =================================
MeshElementTable* bet = model->getMeshBoundaryElements();
outfile << "#group all-boundaries" << endl;
// Loop all boundaries
index = 0;
while (true) {
BodyElement* be = model->getBoundary(index++);
if (be==NULL) break;
int nof_elements = be->getNofMeshElements();
if (nof_elements == 0) {
continue;
}
const beStatus be_status = be->getStatus();
if ( be_status & (BE_DEVIDED | BE_SWAPPED) ) {
continue;
}
//--Boundary name as the group name
outfile << "#group " << be->getName() << endl;
for (i = 0; i < nof_elements; i++) {
int eid = be->getMeshElementId(i);
const int* nodes = bet->getNodeIds(eid);
int parent1_id = bet->parentIds[eid][0];
int body1_id = bt->parentIds[parent1_id][0];
meshElementCode elem_code = bet->getElementCode(eid);
outfile << "all ";
// Elm type for the boundary element
outfile << MeshElementDesc[elem_code][DESC_ELEM_TYPE];
// Internal node ids for the boundary element
int nof_nodes = MeshElementDesc[elem_code][DESC_NOF_NODES];
outfile << " ";
for (int k = 0; k < nof_nodes; k++) {
outfile << " " << nodes[k];
}
outfile << endl;
}
//--End for this boundary
outfile << "#endgroup " << be->getName() << endl;
}
//---End for all boundaries
outfile << "#endgroup all-boundaries" << endl;
//---End for all elements
outfile << "#endgroup all" << endl;
// All done
// ========
gui->showMsg("Mesh saved in Elmer Post format!");
}
// Thetis format file
// ==================
void
ModelOutputManager::write_Thetis_mesh(ostream& outfile)
{
const ModelInfo* mi = model->getModelInfo();
const MeshInfo* mei = model->getMeshInfo();
const ModelStatistics* ms = model->getModelStatistics();
UserInterface* gui = theControlCenter->getGui();
int i,j, index;
int nofNodes = mei->nofNodes;
int nofBulkElements = mei->nofBulkElements;
int nofBulkElements_net = mei->nofBulkElements - mei->nofSplittedElements;
int nofBodies = mei->nofBodies;
int nofBoundaries = ms->nofInnerBoundaries + ms->nofOuterBoundaries;
int nofBndrElements = mei->nofBoundaryElements;
//---Header
outfile << "# nof: Nodes, Max external-id, Elements, Max external-id, Bodies, Nof boundaries, Nof boundary elements\n";
outfile << nofNodes << " " << mei->maxExternalNodeId;
outfile << " " << nofBulkElements_net << " " << mei->maxExternalElementId;
outfile << " " << nofBodies;
outfile << " " << nofBoundaries;
outfile << " " << nofBndrElements;
outfile << endl;
//---Node coordinates
outfile << "##" << endl;
outfile << "#Nodes:" << endl;
//outfile << "External-id, Coordinates, Nof bodies, External body ids" << endl;
outfile << "# External-id, Coordinates" << endl;
Point3* nodes = model->getMeshNodeData();
for (i = 0; i < nofNodes; i++) {
Point3& p = nodes[i];
outfile << model->getMeshNodeIdInt2Ext(i); // External id
outfile << " " << p[0] << " " << p[1] << " " << p[2];
// nof node bodies and external body_ids
//short nof_node_bodies = meshNode2Bodies[i][0];
//outfile << " " << nof_node_bodies;
//for (short j = 0; j < nof_node_bodies; j++) {
//int body_id = meshNode2Bodies[i][1 + j];
//outfile << " " << meshBodyInt2Ext[body_id];
//}
outfile << endl;
}
MeshElementTable* bt = model->getMeshBulkElements();
MeshElementTable* bet = model->getMeshBoundaryElements();
//---Elements
outfile << "##" << endl;
outfile << "#Elements:" << endl;
outfile << "# External-id, External body-id, Element-code" << endl;
outfile << "# Nof nodes, External node-ids" << endl;
outfile << "# Nof neighbour elements , External neighbour ids" << endl;
for (i = 0; i < nofBulkElements; i++) {
if ( bt->splitted[i] ) {
continue;
}
const int* node_ids = bt->getNodeIds(i);
meshElementCode elem_code = bt->getElementCode(i);
int body_id = bt->parentIds[i][0];
int body_tag = model->getModelObjectTagById(body_id);
outfile << model->getMeshBulkElementIdInt2Ext(i);
outfile << " " << model->getBodyTagInt2Ext(body_tag);
outfile << " " << MeshElementDesc[elem_code][DESC_ELEM_TYPE];
// External node ids
outfile << " ";
int nof_nodes = MeshElementDesc[elem_code][DESC_NOF_NODES];
outfile << nof_nodes << " ";
for (j = 0; j < nof_nodes; j++) {
outfile << " " << model->getMeshNodeIdInt2Ext(node_ids[j]);
}
// External neigbour element ids
outfile << " ";
int nof_bndr_elems = MeshElementDesc[elem_code][DESC_NOF_BNDR_ELEMS];
outfile << nof_bndr_elems << " ";
for (j = 0; j < nof_bndr_elems; j++) {
int elem_id = bt->neighborIds[i][j];
outfile << " " << model->getMeshBulkElementIdInt2Ext(elem_id);
}
outfile << endl;
}
//---Boundaries
outfile << "##" << endl;
outfile << "#Boundaries:" << endl;
outfile << "# Boundary id, "
<< "External body1-id, body2-id, Nof bndr elements" << endl;
index = 0;
while (true) {
BodyElement* be = model->getBoundary(index++);
if (be==NULL) break;
// Name as comment line
outfile << "#" << be->getName() << endl;
// Parent tags
int body1_tag = be->getParentTag(1);
int body2_tag = be->getParentTag(2);
outfile << be->Tag();
outfile << " " << model->getBodyTagInt2Ext(body1_tag);
if (body2_tag != NO_INDEX)
outfile << " " << model->getBodyTagInt2Ext(body2_tag);
else
outfile << " " << NO_INDEX;
// Nof elements
outfile << " " << be->getNofMeshElements();
outfile << endl;
}
//---Boundary elements
outfile << "##" << endl;
outfile << "#Boundary elements:" << endl;
outfile << "# Sequence nbr, Boundary id" << endl;
outfile << "# External parent-elem1-id, parent-elem2-id" << endl;
outfile << "# External body1-id, body2-id, Element-code" << endl;
outfile << "# Nof nodes, External node-ids" << endl;
// Loop all boundary elements
index = 0;
while (true) {
BodyElement* be = model->getBoundary(index++);
if (be==NULL) break;
int nof_elements = be->getNofMeshElements();
for (i = 0; i < nof_elements; i++) {
int index = be->getMeshElementId(i);
int bndr_tag = be->Tag();
const int* node_ids = bet->getNodeIds(index);
meshElementCode elem_code = bet->getElementCode(index);
outfile << i + 1 // Sequence nbr
<< " "
<< bndr_tag; // Boundary tag
int p1_id = bet->parentIds[index][0]; // Int parent elem1
int p2_id = bet->parentIds[index][1]; // Int parent elem2
// parent elements
outfile << " " << model->getMeshBulkElementIdInt2Ext(p1_id); // -->ext parent elem1
outfile << " " << model->getMeshBulkElementIdInt2Ext(p2_id); // -->ext parent elem2
// parent bodies
int b1_id = bt->parentIds[p1_id][0]; // Int parent body1
int b1_tag = model->getModelObjectTagById(b1_id);
outfile << " " << model->getBodyTagInt2Ext(b1_tag); // Ext parent body1
if (p2_id != NO_INDEX) {
int b2_id = bt->parentIds[p2_id][0];// Int parent body2
int b2_tag = model->getModelObjectTagById(b2_id);
outfile << " " << model->getBodyTagInt2Ext(b2_tag); // Ext parent body2
} else {
outfile << " " << NO_INDEX; // No parent elem2
}
// Elm type for the boundary element
outfile << " " << MeshElementDesc[elem_code][DESC_ELEM_TYPE];
// Nodes in the boundary element
int nof_nodes = MeshElementDesc[elem_code][DESC_NOF_NODES];
outfile << " ";
outfile << nof_nodes << " ";
// External node ids for the boundary element
for (int k = 0; k < nof_nodes; k++) {
outfile << " " << model->getMeshNodeIdInt2Ext(node_ids[k]);
}
outfile << endl;
}
}
gui->showMsg("Mesh saved in Thetis format!");
}
syntax highlighted by Code2HTML, v. 0.9.1