weipa/src/FinleyDomain.cpp


/*******************************************************
*
* Copyright (c) 2003-2010 by University of Queensland
* Earth Systems Science Computational Center (ESSCC)
* http://www.uq.edu.au/esscc
*
* Primary Business: Queensland, Australia
* Licensed under the Open Software License version 3.0
* http://www.opensource.org/licenses/osl-3.0.php
*
*******************************************************/

#include <weipa/FinleyDomain.h>
#include <weipa/FinleyNodes.h>
#include <weipa/DataVar.h>

#ifndef VISIT_PLUGIN
#include <dudley/CppAdapter/MeshAdapter.h>
#include <finley/CppAdapter/MeshAdapter.h>
extern "C" {
#include <dudley/Mesh.h>
#include <finley/Mesh.h>
}
#endif

#include <iostream>

#if USE_NETCDF
#include <netcdfcpp.h>
#endif

#if USE_SILO
#include <silo.h>
#endif

using namespace std;

namespace weipa {

//
//
//
FinleyDomain::FinleyDomain() :
    initialized(false)
{
}

//
//
//
FinleyDomain::FinleyDomain(const FinleyDomain& m) :
    boost::enable_shared_from_this<FinleyDomain>()
{
    nodes = FinleyNodes_ptr(new FinleyNodes(*m.nodes));
    cells = FinleyElements_ptr(new FinleyElements(*m.cells));
    faces = FinleyElements_ptr(new FinleyElements(*m.faces));
    contacts = FinleyElements_ptr(new FinleyElements(*m.contacts));
    initialized = m.initialized;
}

//
//
//
FinleyDomain::~FinleyDomain()
{
    cleanup();
}

//
//
//
void FinleyDomain::cleanup()
{
    nodes.reset();
    cells.reset();
    faces.reset();
    contacts.reset();
    initialized = false;
}

//
//
//
bool FinleyDomain::initFromEscript(const escript::AbstractDomain* escriptDomain)
{
#ifndef VISIT_PLUGIN
    cleanup();

    if (dynamic_cast<const finley::MeshAdapter*>(escriptDomain)) {
        const Finley_Mesh* finleyMesh =
            dynamic_cast<const finley::MeshAdapter*>(escriptDomain)
                ->getFinley_Mesh();

        nodes = FinleyNodes_ptr(new FinleyNodes("Elements"));
        cells = FinleyElements_ptr(new FinleyElements("Elements", nodes));
        faces = FinleyElements_ptr(new FinleyElements("FaceElements", nodes));
        contacts = FinleyElements_ptr(new FinleyElements("ContactElements", nodes));

        if (nodes->initFromFinley(finleyMesh->Nodes) &&
                cells->initFromFinley(finleyMesh->Elements) &&
                faces->initFromFinley(finleyMesh->FaceElements) &&
                contacts->initFromFinley(finleyMesh->ContactElements)) {
            initialized = true;
        }
    } else if (dynamic_cast<const dudley::MeshAdapter*>(escriptDomain)) {
        const Dudley_Mesh* dudleyMesh =
            dynamic_cast<const dudley::MeshAdapter*>(escriptDomain)
                ->getDudley_Mesh();

        nodes = FinleyNodes_ptr(new FinleyNodes("Elements"));
        cells = FinleyElements_ptr(new FinleyElements("Elements", nodes));
        faces = FinleyElements_ptr(new FinleyElements("FaceElements", nodes));
        contacts = FinleyElements_ptr(new FinleyElements("ContactElements", nodes));

        if (nodes->initFromDudley(dudleyMesh->Nodes) &&
                cells->initFromDudley(dudleyMesh->Elements) &&
                faces->initFromDudley(dudleyMesh->FaceElements)) {
            initialized = true;
        }
    }

    return initialized;
#else // VISIT_PLUGIN
    return false;
#endif
}

//
// Reads mesh and element data from NetCDF file with given name
//
bool FinleyDomain::initFromFile(const string& filename)
{
    cleanup();
    
#if USE_NETCDF
    NcError ncerr(NcError::silent_nonfatal);
    NcFile* input;
 
    input = new NcFile(filename.c_str());
    if (!input->is_valid()) {
        cerr << "Could not open input file " << filename << "." << endl;
        delete input;
        return false;
    }

    nodes = FinleyNodes_ptr(new FinleyNodes("Elements"));
    if (!nodes->readFromNc(input))
        return false;

    // Read all element types
    cells = FinleyElements_ptr(new FinleyElements("Elements", nodes));
    cells->readFromNc(input);
    faces = FinleyElements_ptr(new FinleyElements("FaceElements", nodes));
    faces->readFromNc(input);
    contacts = FinleyElements_ptr(new FinleyElements("ContactElements", nodes));
    contacts->readFromNc(input);

    delete input;
    initialized = true;
#endif

    return initialized;
}

Centering FinleyDomain::getCenteringForFunctionSpace(int fsCode) const
{
    return (fsCode==FINLEY_REDUCED_NODES || fsCode==FINLEY_NODES ?
        NODE_CENTERED : ZONE_CENTERED);
}

//
//
//
NodeData_ptr FinleyDomain::getMeshForFunctionSpace(int fsCode) const
{
    NodeData_ptr result;

    if (!initialized)
        return result;

    ElementData_ptr elements = getElementsForFunctionSpace(fsCode);
    if (elements != NULL)
        result = elements->getNodes();
 
    return result;
}

//
//
//
ElementData_ptr FinleyDomain::getElementsForFunctionSpace(int fsCode) const
{
    ElementData_ptr result;

    if (!initialized) {
        return result;
    }

    if (fsCode == FINLEY_NODES) {
        result = cells;
    } else if (fsCode == FINLEY_REDUCED_NODES) {
        result = cells->getReducedElements();
        if (!result)
            result = cells;
    } else {
        switch (fsCode) {
            case FINLEY_REDUCED_ELEMENTS:
            case FINLEY_ELEMENTS:
                result = cells;
                break;

            case FINLEY_REDUCED_FACE_ELEMENTS:
            case FINLEY_FACE_ELEMENTS:
                result = faces;
                break;

            case FINLEY_REDUCED_CONTACT_ELEMENTS_1:
            case FINLEY_REDUCED_CONTACT_ELEMENTS_2:
            case FINLEY_CONTACT_ELEMENTS_1:
            case FINLEY_CONTACT_ELEMENTS_2:
                result = contacts;
                break;

            default: {
                cerr << "Unsupported function space type " << fsCode
                    << "!" << endl;
                return result;
            }
        }
        int typeId = static_cast<FinleyElements*>(result.get())
            ->getFinleyTypeId();
        if (typeId != Finley_Line3Macro &&
                typeId != Finley_Rec9 && typeId != Finley_Rec9Macro &&
                typeId != Finley_Hex27 &&
                typeId != Finley_Hex27Macro && typeId != Finley_Tri6 &&
                typeId != Finley_Tri6Macro && typeId != Finley_Tet10 &&
                typeId != Finley_Tet10Macro) {
            if (result->getReducedElements())
                result = result->getReducedElements();
        }
    }

    return result;
}

//
// Returns a vector of strings containing mesh names for this domain
//
StringVec FinleyDomain::getMeshNames() const
{
    StringVec res;
    if (initialized) {
        StringVec tmpVec;
        tmpVec = cells->getMeshNames();
        res.insert(res.end(), tmpVec.begin(), tmpVec.end());
        tmpVec = faces->getMeshNames();
        res.insert(res.end(), tmpVec.begin(), tmpVec.end());
        tmpVec = contacts->getMeshNames();
        res.insert(res.end(), tmpVec.begin(), tmpVec.end());
    }
    return res;
}

//
// Returns a vector of strings containing mesh variable names for this domain
//
StringVec FinleyDomain::getVarNames() const
{
    StringVec res;
 
    if (initialized) {
        res = nodes->getVarNames();
        StringVec tmpVec = cells->getVarNames();
        res.insert(res.end(), tmpVec.begin(), tmpVec.end());
        tmpVec = faces->getVarNames();
        res.insert(res.end(), tmpVec.begin(), tmpVec.end());
        tmpVec = contacts->getVarNames();
        res.insert(res.end(), tmpVec.begin(), tmpVec.end());
    }

    return res;
}

//
//
//
DataVar_ptr FinleyDomain::getDataVarByName(const string& name) const
{
    if (!initialized) {
        throw "Domain not initialized";
    }

    DataVar_ptr var(new DataVar(name));
    if (name.find("ContactElements_") != name.npos) {
        const IntVec& data = contacts->getVarDataByName(name);
        string elementName = name.substr(0, name.find('_'));
        ElementData_ptr elements = getElementsByName(elementName);
        var->initFromMeshData(shared_from_this(), data,
                FINLEY_CONTACT_ELEMENTS_1, ZONE_CENTERED, elements->getNodes(),
                elements->getIDs());
    } else if (name.find("FaceElements_") != name.npos) {
        const IntVec& data =  faces->getVarDataByName(name);
        string elementName = name.substr(0, name.find('_'));
        ElementData_ptr elements = getElementsByName(elementName);
        var->initFromMeshData(shared_from_this(), data,
                FINLEY_FACE_ELEMENTS, ZONE_CENTERED, elements->getNodes(),
                elements->getIDs());
    } else if (name.find("Elements_") != name.npos) {
        const IntVec& data =  cells->getVarDataByName(name);
        string elementName = name.substr(0, name.find('_'));
        ElementData_ptr elements = getElementsByName(elementName);
        var->initFromMeshData(shared_from_this(), data, FINLEY_ELEMENTS,
                ZONE_CENTERED, elements->getNodes(), elements->getIDs());
    } else if (name.find("Nodes_") != name.npos) {
        const IntVec& data =  nodes->getVarDataByName(name);
        var->initFromMeshData(shared_from_this(), data, FINLEY_NODES,
                NODE_CENTERED, getNodes(), getNodes()->getNodeIDs());
    } else {
        cerr << "WARNING: Unrecognized domain variable '" << name << "'\n";
        return DataVar_ptr();
    }

    return var;
}

//
//
//
ElementData_ptr FinleyDomain::getElementsByName(const string& name) const
{
    ElementData_ptr ret;
    if (name == "Elements")
        ret = cells;
    else if (name == "ReducedElements")
        ret = cells->getReducedElements();
    else if (name == "FaceElements")
        ret = faces;
    else if (name == "ReducedFaceElements")
        ret = faces->getReducedElements();
    else if (name == "ContactElements")
        ret = contacts;
    else if (name == "ReducedContactElements")
        ret = contacts->getReducedElements();

    return ret;
}

//
//
//
NodeData_ptr FinleyDomain::getMeshByName(const string& name) const
{
    NodeData_ptr ret;
    if (initialized) {
        ElementData_ptr els = getElementsByName(name);
        if (els)
            ret = els->getNodes();
    }

    return ret;
}

//
//
//
void FinleyDomain::reorderGhostZones(int ownIndex)
{
    if (initialized) {
        cells->reorderGhostZones(ownIndex);
        faces->reorderGhostZones(ownIndex);
        contacts->reorderGhostZones(ownIndex);
#ifdef _DEBUG
        cout << "block " << ownIndex << " has " << cells->getGhostCount()
             << " ghost zones," << endl;
        cout << "\t" << faces->getGhostCount() << " ghost faces," << endl;
        cout << "\t" << contacts->getGhostCount() << " ghost contacts." << endl;
#endif
    }
}

//
//
//
void FinleyDomain::removeGhostZones(int ownIndex)
{
    if (initialized) {
        cells->removeGhostZones(ownIndex);
        faces->removeGhostZones(ownIndex);
        contacts->removeGhostZones(ownIndex);
#ifdef _DEBUG
        cout << "After removing ghost zones there are" << endl;
        cout << "    " << nodes->getNumNodes() << " Nodes, ";
        cout << cells->getCount() << " Elements, ";
        cout << faces->getCount() << " Face elements, ";
        cout << contacts->getCount() << " Contact elements left." << endl;
#endif
    }
}

//
//
//
bool FinleyDomain::writeToSilo(DBfile* dbfile, const string& pathInSilo,
                               const StringVec& labels, const StringVec& units)
{
#if USE_SILO
    // Write nodes, elements and mesh variables
    if (!initialized ||
            !cells->writeToSilo(dbfile, pathInSilo, labels, units) ||
            !faces->writeToSilo(dbfile, pathInSilo, labels, units) ||
            !contacts->writeToSilo(dbfile, pathInSilo, labels, units))
        return false;

    siloPath = pathInSilo;
    return true;

#else // !USE_SILO
    return false;
#endif
}

} // namespace weipa

1
2	/*******************************************************
3	*
4	* Copyright (c) 2003-2010 by University of Queensland
5	* Earth Systems Science Computational Center (ESSCC)
6	* http://www.uq.edu.au/esscc
7	*
8	* Primary Business: Queensland, Australia
9	* Licensed under the Open Software License version 3.0
10	* http://www.opensource.org/licenses/osl-3.0.php
11	*
12	*******************************************************/
13
14	#include <weipa/FinleyDomain.h>
15	#include <weipa/FinleyNodes.h>
16	#include <weipa/DataVar.h>
17
18	#ifndef VISIT_PLUGIN
19	#include <dudley/CppAdapter/MeshAdapter.h>
20	#include <finley/CppAdapter/MeshAdapter.h>
21	extern "C" {
22	#include <dudley/Mesh.h>
23	#include <finley/Mesh.h>
24	}
25	#endif
26
27	#include <iostream>
28
29	#if USE_NETCDF
30	#include <netcdfcpp.h>
31	#endif
32
33	#if USE_SILO
34	#include <silo.h>
35	#endif
36
37	using namespace std;
38
39	namespace weipa {
40
41	//
42	//
43	//
44	FinleyDomain::FinleyDomain() :
45	initialized(false)
46	{
47	}
48
49	//
50	//
51	//
52	FinleyDomain::FinleyDomain(const FinleyDomain& m) :
53	boost::enable_shared_from_this<FinleyDomain>()
54	{
55	nodes = FinleyNodes_ptr(new FinleyNodes(*m.nodes));
56	cells = FinleyElements_ptr(new FinleyElements(*m.cells));
57	faces = FinleyElements_ptr(new FinleyElements(*m.faces));
58	contacts = FinleyElements_ptr(new FinleyElements(*m.contacts));
59	initialized = m.initialized;
60	}
61
62	//
63	//
64	//
65	FinleyDomain::~FinleyDomain()
66	{
67	cleanup();
68	}
69
70	//
71	//
72	//
73	void FinleyDomain::cleanup()
74	{
75	nodes.reset();
76	cells.reset();
77	faces.reset();
78	contacts.reset();
79	initialized = false;
80	}
81
82	//
83	//
84	//
85	bool FinleyDomain::initFromEscript(const escript::AbstractDomain* escriptDomain)
86	{
87	#ifndef VISIT_PLUGIN
88	cleanup();
89
90	if (dynamic_cast<const finley::MeshAdapter*>(escriptDomain)) {
91	const Finley_Mesh* finleyMesh =
92	dynamic_cast<const finley::MeshAdapter*>(escriptDomain)
93	->getFinley_Mesh();
94
95	nodes = FinleyNodes_ptr(new FinleyNodes("Elements"));
96	cells = FinleyElements_ptr(new FinleyElements("Elements", nodes));
97	faces = FinleyElements_ptr(new FinleyElements("FaceElements", nodes));
98	contacts = FinleyElements_ptr(new FinleyElements("ContactElements", nodes));
99
100	if (nodes->initFromFinley(finleyMesh->Nodes) &&
101	cells->initFromFinley(finleyMesh->Elements) &&
102	faces->initFromFinley(finleyMesh->FaceElements) &&
103	contacts->initFromFinley(finleyMesh->ContactElements)) {
104	initialized = true;
105	}
106	} else if (dynamic_cast<const dudley::MeshAdapter*>(escriptDomain)) {
107	const Dudley_Mesh* dudleyMesh =
108	dynamic_cast<const dudley::MeshAdapter*>(escriptDomain)
109	->getDudley_Mesh();
110
111	nodes = FinleyNodes_ptr(new FinleyNodes("Elements"));
112	cells = FinleyElements_ptr(new FinleyElements("Elements", nodes));
113	faces = FinleyElements_ptr(new FinleyElements("FaceElements", nodes));
114	contacts = FinleyElements_ptr(new FinleyElements("ContactElements", nodes));
115
116	if (nodes->initFromDudley(dudleyMesh->Nodes) &&
117	cells->initFromDudley(dudleyMesh->Elements) &&
118	faces->initFromDudley(dudleyMesh->FaceElements)) {
119	initialized = true;
120	}
121	}
122
123	return initialized;
124	#else // VISIT_PLUGIN
125	return false;
126	#endif
127	}
128
129	//
130	// Reads mesh and element data from NetCDF file with given name
131	//
132	bool FinleyDomain::initFromFile(const string& filename)
133	{
134	cleanup();
135
136	#if USE_NETCDF
137	NcError ncerr(NcError::silent_nonfatal);
138	NcFile* input;
139
140	input = new NcFile(filename.c_str());
141	if (!input->is_valid()) {
142	cerr << "Could not open input file " << filename << "." << endl;
143	delete input;
144	return false;
145	}
146
147	nodes = FinleyNodes_ptr(new FinleyNodes("Elements"));
148	if (!nodes->readFromNc(input))
149	return false;
150
151	// Read all element types
152	cells = FinleyElements_ptr(new FinleyElements("Elements", nodes));
153	cells->readFromNc(input);
154	faces = FinleyElements_ptr(new FinleyElements("FaceElements", nodes));
155	faces->readFromNc(input);
156	contacts = FinleyElements_ptr(new FinleyElements("ContactElements", nodes));
157	contacts->readFromNc(input);
158
159	delete input;
160	initialized = true;
161	#endif
162
163	return initialized;
164	}
165
166	Centering FinleyDomain::getCenteringForFunctionSpace(int fsCode) const
167	{
168	return (fsCode==FINLEY_REDUCED_NODES \|\| fsCode==FINLEY_NODES ?
169	NODE_CENTERED : ZONE_CENTERED);
170	}
171
172	//
173	//
174	//
175	NodeData_ptr FinleyDomain::getMeshForFunctionSpace(int fsCode) const
176	{
177	NodeData_ptr result;
178
179	if (!initialized)
180	return result;
181
182	ElementData_ptr elements = getElementsForFunctionSpace(fsCode);
183	if (elements != NULL)
184	result = elements->getNodes();
185
186	return result;
187	}
188
189	//
190	//
191	//
192	ElementData_ptr FinleyDomain::getElementsForFunctionSpace(int fsCode) const
193	{
194	ElementData_ptr result;
195
196	if (!initialized) {
197	return result;
198	}
199
200	if (fsCode == FINLEY_NODES) {
201	result = cells;
202	} else if (fsCode == FINLEY_REDUCED_NODES) {
203	result = cells->getReducedElements();
204	if (!result)
205	result = cells;
206	} else {
207	switch (fsCode) {
208	case FINLEY_REDUCED_ELEMENTS:
209	case FINLEY_ELEMENTS:
210	result = cells;
211	break;
212
213	case FINLEY_REDUCED_FACE_ELEMENTS:
214	case FINLEY_FACE_ELEMENTS:
215	result = faces;
216	break;
217
218	case FINLEY_REDUCED_CONTACT_ELEMENTS_1:
219	case FINLEY_REDUCED_CONTACT_ELEMENTS_2:
220	case FINLEY_CONTACT_ELEMENTS_1:
221	case FINLEY_CONTACT_ELEMENTS_2:
222	result = contacts;
223	break;
224
225	default: {
226	cerr << "Unsupported function space type " << fsCode
227	<< "!" << endl;
228	return result;
229	}
230	}
231	int typeId = static_cast<FinleyElements*>(result.get())
232	->getFinleyTypeId();
233	if (typeId != Finley_Line3Macro &&
234	typeId != Finley_Rec9 && typeId != Finley_Rec9Macro &&
235	typeId != Finley_Hex27 &&
236	typeId != Finley_Hex27Macro && typeId != Finley_Tri6 &&
237	typeId != Finley_Tri6Macro && typeId != Finley_Tet10 &&
238	typeId != Finley_Tet10Macro) {
239	if (result->getReducedElements())
240	result = result->getReducedElements();
241	}
242	}
243
244	return result;
245	}
246
247	//
248	// Returns a vector of strings containing mesh names for this domain
249	//
250	StringVec FinleyDomain::getMeshNames() const
251	{
252	StringVec res;
253	if (initialized) {
254	StringVec tmpVec;
255	tmpVec = cells->getMeshNames();
256	res.insert(res.end(), tmpVec.begin(), tmpVec.end());
257	tmpVec = faces->getMeshNames();
258	res.insert(res.end(), tmpVec.begin(), tmpVec.end());
259	tmpVec = contacts->getMeshNames();
260	res.insert(res.end(), tmpVec.begin(), tmpVec.end());
261	}
262	return res;
263	}
264
265	//
266	// Returns a vector of strings containing mesh variable names for this domain
267	//
268	StringVec FinleyDomain::getVarNames() const
269	{
270	StringVec res;
271
272	if (initialized) {
273	res = nodes->getVarNames();
274	StringVec tmpVec = cells->getVarNames();
275	res.insert(res.end(), tmpVec.begin(), tmpVec.end());
276	tmpVec = faces->getVarNames();
277	res.insert(res.end(), tmpVec.begin(), tmpVec.end());
278	tmpVec = contacts->getVarNames();
279	res.insert(res.end(), tmpVec.begin(), tmpVec.end());
280	}
281
282	return res;
283	}
284
285	//
286	//
287	//
288	DataVar_ptr FinleyDomain::getDataVarByName(const string& name) const
289	{
290	if (!initialized) {
291	throw "Domain not initialized";
292	}
293
294	DataVar_ptr var(new DataVar(name));
295	if (name.find("ContactElements_") != name.npos) {
296	const IntVec& data = contacts->getVarDataByName(name);
297	string elementName = name.substr(0, name.find('_'));
298	ElementData_ptr elements = getElementsByName(elementName);
299	var->initFromMeshData(shared_from_this(), data,
300	FINLEY_CONTACT_ELEMENTS_1, ZONE_CENTERED, elements->getNodes(),
301	elements->getIDs());
302	} else if (name.find("FaceElements_") != name.npos) {
303	const IntVec& data = faces->getVarDataByName(name);
304	string elementName = name.substr(0, name.find('_'));
305	ElementData_ptr elements = getElementsByName(elementName);
306	var->initFromMeshData(shared_from_this(), data,
307	FINLEY_FACE_ELEMENTS, ZONE_CENTERED, elements->getNodes(),
308	elements->getIDs());
309	} else if (name.find("Elements_") != name.npos) {
310	const IntVec& data = cells->getVarDataByName(name);
311	string elementName = name.substr(0, name.find('_'));
312	ElementData_ptr elements = getElementsByName(elementName);
313	var->initFromMeshData(shared_from_this(), data, FINLEY_ELEMENTS,
314	ZONE_CENTERED, elements->getNodes(), elements->getIDs());
315	} else if (name.find("Nodes_") != name.npos) {
316	const IntVec& data = nodes->getVarDataByName(name);
317	var->initFromMeshData(shared_from_this(), data, FINLEY_NODES,
318	NODE_CENTERED, getNodes(), getNodes()->getNodeIDs());
319	} else {
320	cerr << "WARNING: Unrecognized domain variable '" << name << "'\n";
321	return DataVar_ptr();
322	}
323
324	return var;
325	}
326
327	//
328	//
329	//
330	ElementData_ptr FinleyDomain::getElementsByName(const string& name) const
331	{
332	ElementData_ptr ret;
333	if (name == "Elements")
334	ret = cells;
335	else if (name == "ReducedElements")
336	ret = cells->getReducedElements();
337	else if (name == "FaceElements")
338	ret = faces;
339	else if (name == "ReducedFaceElements")
340	ret = faces->getReducedElements();
341	else if (name == "ContactElements")
342	ret = contacts;
343	else if (name == "ReducedContactElements")
344	ret = contacts->getReducedElements();
345
346	return ret;
347	}
348
349	//
350	//
351	//
352	NodeData_ptr FinleyDomain::getMeshByName(const string& name) const
353	{
354	NodeData_ptr ret;
355	if (initialized) {
356	ElementData_ptr els = getElementsByName(name);
357	if (els)
358	ret = els->getNodes();
359	}
360
361	return ret;
362	}
363
364	//
365	//
366	//
367	void FinleyDomain::reorderGhostZones(int ownIndex)
368	{
369	if (initialized) {
370	cells->reorderGhostZones(ownIndex);
371	faces->reorderGhostZones(ownIndex);
372	contacts->reorderGhostZones(ownIndex);
373	#ifdef _DEBUG
374	cout << "block " << ownIndex << " has " << cells->getGhostCount()
375	<< " ghost zones," << endl;
376	cout << "\t" << faces->getGhostCount() << " ghost faces," << endl;
377	cout << "\t" << contacts->getGhostCount() << " ghost contacts." << endl;
378	#endif
379	}
380	}
381
382	//
383	//
384	//
385	void FinleyDomain::removeGhostZones(int ownIndex)
386	{
387	if (initialized) {
388	cells->removeGhostZones(ownIndex);
389	faces->removeGhostZones(ownIndex);
390	contacts->removeGhostZones(ownIndex);
391	#ifdef _DEBUG
392	cout << "After removing ghost zones there are" << endl;
393	cout << " " << nodes->getNumNodes() << " Nodes, ";
394	cout << cells->getCount() << " Elements, ";
395	cout << faces->getCount() << " Face elements, ";
396	cout << contacts->getCount() << " Contact elements left." << endl;
397	#endif
398	}
399	}
400
401	//
402	//
403	//
404	bool FinleyDomain::writeToSilo(DBfile* dbfile, const string& pathInSilo,
405	const StringVec& labels, const StringVec& units)
406	{
407	#if USE_SILO
408	// Write nodes, elements and mesh variables
409	if (!initialized \|\|
410	!cells->writeToSilo(dbfile, pathInSilo, labels, units) \|\|
411	!faces->writeToSilo(dbfile, pathInSilo, labels, units) \|\|
412	!contacts->writeToSilo(dbfile, pathInSilo, labels, units))
413	return false;
414
415	siloPath = pathInSilo;
416	return true;
417
418	#else // !USE_SILO
419	return false;
420	#endif
421	}
422
423	} // namespace weipa
424