From 3185a4c67b74fab036f554d44ed0e5a2f0d69e60 Mon Sep 17 00:00:00 2001 From: schnetter Date: Thu, 8 Apr 2004 16:08:43 +0000 Subject: A symmetry thorn with a 90 degree rotating symmetry. This is not yet tested. git-svn-id: http://svn.cactuscode.org/arrangements/CactusNumerical/RotatingSymmetry90/trunk@2 c3c03602-0f4f-0410-b3fa-d2c81c8a7dc5 --- README | 23 +++ doc/documentation.tex | 143 +++++++++++++++ interface.ccl | 33 ++++ param.ccl | 14 ++ schedule.ccl | 12 ++ src/make.code.defn | 8 + src/registersymmetry.c | 41 +++++ src/rotatingsymmetry90.c | 455 +++++++++++++++++++++++++++++++++++++++++++++++ 8 files changed, 729 insertions(+) create mode 100644 README create mode 100644 doc/documentation.tex create mode 100644 interface.ccl create mode 100644 param.ccl create mode 100644 schedule.ccl create mode 100644 src/make.code.defn create mode 100644 src/registersymmetry.c create mode 100644 src/rotatingsymmetry90.c diff --git a/README b/README new file mode 100644 index 0000000..084d252 --- /dev/null +++ b/README @@ -0,0 +1,23 @@ +CVS info : $Header$ + +Cactus Code Thorn RotatingSymmetry90 +Thorn Author(s) : Erik Schnetter +Thorn Maintainer(s) : Erik Schnetter +-------------------------------------------------------------------------- + +Purpose of the thorn: + +Provide a 90 degree rotational symmetry boundary condition. This is +the kind of symmetry that exists e.g. in a singl black hole evolution +when the black hole has spin. The same setup also admits an +additional bitant symmetry. + +Currently the symmetry is fixed to be located at the lower x and y +boundaries, so that there is a rotation about the z axis. + +The coordinate thorn does not know about this symmetry. You will have +to set up the domain such that you have an upper quadrant in the x and +y directions. + +This symmetry boundary condition can be mixed with arbitrary symmetry +conditions in the z direction, e.g. bitant or periodicity. diff --git a/doc/documentation.tex b/doc/documentation.tex new file mode 100644 index 0000000..23b2c81 --- /dev/null +++ b/doc/documentation.tex @@ -0,0 +1,143 @@ +% *======================================================================* +% Cactus Thorn template for ThornGuide documentation +% Author: Ian Kelley +% Date: Sun Jun 02, 2002 +% $Header$ +% +% Thorn documentation in the latex file doc/documentation.tex +% will be included in ThornGuides built with the Cactus make system. +% The scripts employed by the make system automatically include +% pages about variables, parameters and scheduling parsed from the +% relevant thorn CCL files. +% +% This template contains guidelines which help to assure that your +% documentation will be correctly added to ThornGuides. More +% information is available in the Cactus UsersGuide. +% +% Guidelines: +% - Do not change anything before the line +% % START CACTUS THORNGUIDE", +% except for filling in the title, author, date, etc. fields. +% - Each of these fields should only be on ONE line. +% - Author names should be separated with a \\ or a comma. +% - You can define your own macros, but they must appear after +% the START CACTUS THORNGUIDE line, and must not redefine standard +% latex commands. +% - To avoid name clashes with other thorns, 'labels', 'citations', +% 'references', and 'image' names should conform to the following +% convention: +% ARRANGEMENT_THORN_LABEL +% For example, an image wave.eps in the arrangement CactusWave and +% thorn WaveToyC should be renamed to CactusWave_WaveToyC_wave.eps +% - Graphics should only be included using the graphicx package. +% More specifically, with the "\includegraphics" command. Do +% not specify any graphic file extensions in your .tex file. This +% will allow us to create a PDF version of the ThornGuide +% via pdflatex. +% - References should be included with the latex "\bibitem" command. +% - Use \begin{abstract}...\end{abstract} instead of \abstract{...} +% - Do not use \appendix, instead include any appendices you need as +% standard sections. +% - For the benefit of our Perl scripts, and for future extensions, +% please use simple latex. +% +% *======================================================================* +% +% Example of including a graphic image: +% \begin{figure}[ht] +% \begin{center} +% \includegraphics[width=6cm]{MyArrangement_MyThorn_MyFigure} +% \end{center} +% \caption{Illustration of this and that} +% \label{MyArrangement_MyThorn_MyLabel} +% \end{figure} +% +% Example of using a label: +% \label{MyArrangement_MyThorn_MyLabel} +% +% Example of a citation: +% \cite{MyArrangement_MyThorn_Author99} +% +% Example of including a reference +% \bibitem{MyArrangement_MyThorn_Author99} +% {J. Author, {\em The Title of the Book, Journal, or periodical}, 1 (1999), +% 1--16. {\tt http://www.nowhere.com/}} +% +% *======================================================================* + +% If you are using CVS use this line to give version information +% $Header$ + +\documentclass{article} + +% Use the Cactus ThornGuide style file +% (Automatically used from Cactus distribution, if you have a +% thorn without the Cactus Flesh download this from the Cactus +% homepage at www.cactuscode.org) +\usepackage{../../../../doc/latex/cactus} + +\begin{document} + +% The author of the documentation +\author{Erik Schnetter \textless schnetter@aei.mpg.de\textgreater} + +% The title of the document (not necessarily the name of the Thorn) +\title{RotatingSymmetry180} + +% the date your document was last changed, if your document is in CVS, +% please use: +\date{$ $Date$ $} + +\maketitle + +% Do not delete next line +% START CACTUS THORNGUIDE + +% Add all definitions used in this documentation here +% \def\mydef etc + +% Add an abstract for this thorn's documentation +\begin{abstract} + +\end{abstract} + +% The following sections are suggestive only. +% Remove them or add your own. + +\section{Introduction} + +\section{Physical System} + +\section{Numerical Implementation} + +\section{Using This Thorn} + +\subsection{Obtaining This Thorn} + +\subsection{Basic Usage} + +\subsection{Special Behaviour} + +\subsection{Interaction With Other Thorns} + +\subsection{Examples} + +\subsection{Support and Feedback} + +\section{History} + +\subsection{Thorn Source Code} + +\subsection{Thorn Documentation} + +\subsection{Acknowledgements} + + +\begin{thebibliography}{9} + +\end{thebibliography} + +% Do not delete next line +% END CACTUS THORNGUIDE + +\end{document} diff --git a/interface.ccl b/interface.ccl new file mode 100644 index 0000000..46383d9 --- /dev/null +++ b/interface.ccl @@ -0,0 +1,33 @@ +# Interface definition for thorn RotatingSymmetry90 +# $Header$ + +IMPLEMENTS: RotatingSymmetry90 + +INHERITS: Slab + +USES INCLUDE HEADER: Slab.h + + + +CCTK_INT FUNCTION SymmetryRegister (CCTK_STRING IN sym_name) +USES FUNCTION SymmetryRegister + +CCTK_INT FUNCTION \ + SymmetryRegisterGrid \ + (CCTK_POINTER IN cctkGH, \ + CCTK_INT IN sym_handle, \ + CCTK_INT IN ARRAY which_faces, \ + CCTK_INT IN ARRAY symmetry_zone_width) +USES FUNCTION SymmetryRegisterGrid + + + +CCTK_INT FUNCTION Boundary_SelectedGVs \ + (CCTK_POINTER_TO_CONST IN GH, \ + CCTK_INT IN array_size, \ + CCTK_INT ARRAY OUT var_indicies, \ + CCTK_INT ARRAY OUT faces, \ + CCTK_INT ARRAY OUT boundary_widths, \ + CCTK_INT ARRAY OUT table_handles, \ + CCTK_STRING IN bc_name) +USES FUNCTION Boundary_SelectedGVs diff --git a/param.ccl b/param.ccl new file mode 100644 index 0000000..544d191 --- /dev/null +++ b/param.ccl @@ -0,0 +1,14 @@ +# Parameter definitions for thorn RotatingSymmetry90 +# $Header$ + +BOOLEAN verbose "Produce screen output while applying boundary conditions" +{ +} "no" + +BOOLEAN rotating_symmetry_90 "Apply 90 degree rotating symmetry boundary conditions?" +{ +} "no" + +BOOLEAN avoid_origin "Avoid the origin?" +{ +} "no" diff --git a/schedule.ccl b/schedule.ccl new file mode 100644 index 0000000..d624033 --- /dev/null +++ b/schedule.ccl @@ -0,0 +1,12 @@ +# Schedule definitions for thorn RotatingSymmetry90 +# $Header$ + +SCHEDULE Rot90_RegisterSymmetry IN SymmetryRegister +{ + LANG: C +} "Register symmetry boundaries" + +SCHEDULE Rot90_ApplyBC IN BoundaryConditions +{ + LANG: C +} "Apply 90 degree rotational symmetry boundary condition" diff --git a/src/make.code.defn b/src/make.code.defn new file mode 100644 index 0000000..4b71cce --- /dev/null +++ b/src/make.code.defn @@ -0,0 +1,8 @@ +# Main make.code.defn file for thorn RotatingSymmetry90 +# $Header$ + +# Source files in this directory +SRCS = rotatingsymmetry90.c registersymmetry.c + +# Subdirectories containing source files +SUBDIRS = diff --git a/src/registersymmetry.c b/src/registersymmetry.c new file mode 100644 index 0000000..be8999c --- /dev/null +++ b/src/registersymmetry.c @@ -0,0 +1,41 @@ +/* $Header$ */ + +#include "cctk.h" +#include "cctk_Arguments.h" +#include "cctk_Parameters.h" + +void Rot90_RegisterSymmetry (CCTK_ARGUMENTS) +{ + DECLARE_CCTK_ARGUMENTS; + DECLARE_CCTK_PARAMETERS; + + int f; + CCTK_INT handle; + CCTK_INT faces[6]; + CCTK_INT width[6]; + CCTK_INT ierr; + + if (rotating_symmetry_90) { + + for (f=0; f<6; ++f) { + faces[f] = 0; + width[f] = 0; + } + + faces[0] = 1; + width[0] = cctk_nghostzones[0]; + faces[1] = 1; + width[1] = cctk_nghostzones[1]; + + handle = SymmetryRegister ("rotating_symmetry_90"); + if (handle < 0) { + CCTK_WARN (0, "Could not register symmetry boundary condition"); + } + + ierr = SymmetryRegisterGrid (cctkGH, handle, faces, width); + if (ierr < 0) { + CCTK_WARN (0, "Could not register the symmetry boundaries -- probably some other thorn has already registered the same boundary faces for a different symmetry"); + } + + } +} diff --git a/src/rotatingsymmetry90.c b/src/rotatingsymmetry90.c new file mode 100644 index 0000000..f31bb66 --- /dev/null +++ b/src/rotatingsymmetry90.c @@ -0,0 +1,455 @@ +/* $Header$ */ + +#include +#include +#include +#include + +#include "cctk.h" +#include "cctk_Arguments.h" +#include "cctk_Parameters.h" + +#include "util_ErrorCodes.h" +#include "util_Table.h" + +#include "Slab.h" + + + +int BndRot90VI (cGH const * restrict const cctkGH, + int const * restrict const stencil, + int const vi) +{ + DECLARE_CCTK_ARGUMENTS; + DECLARE_CCTK_PARAMETERS; + + int gi; + cGroup group; + cGroupDynamicData data; + char * restrict fullname; + void * restrict varptr; + int firstvar; + + char tensortypealias[100]; + int numvars; + int index; + int srcvi; + void const * restrict srcptr; + int parity; + + int global_bbox[6]; + int global_lbnd[3], global_ubnd[3]; + int fake_bbox[6]; + + CCTK_REAL x0[3], dx[3]; + CCTK_REAL origin[3], dorigin[3]; + int avoid_origin[3]; + int offset[3]; /* offset 0..1 due to avoid_origin */ + + struct xferinfo * restrict xferinfo; + + int have_global_bbox, have_local_bbox; + + int alldirs[2]; + int step; + int ndirs; /* 1 for face, 2 for edge */ + int dir[3]; /* current face or edge */ + int otherdir[3]; /* the other direction of the rotation */ + int q; /* 0..ndirs-1 */ + int d; /* 0..group.dim-1 */ + int ierr; + + /* Check arguments */ + assert (cctkGH); + assert (stencil); + assert (vi>=0 && vi=0 && gi0); + firstvar = CCTK_FirstVarIndexI(gi); + assert (firstvar>=0); + index = vi - firstvar; + assert (index>=0 && index=0); + + ierr = Util_TableGetString + (table, sizeof tensortypealias, tensortypealias, "tensortypealias"); + if (ierr == UTIL_ERROR_TABLE_NO_SUCH_KEY) { + char * groupname = CCTK_GroupName(gi); + assert (groupname); + CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, + "Tensor type alias not declared for group \"%s\" -- assuming a scalar", + groupname); + free (groupname); + strcpy (tensortypealias, "scalar"); + } else if (ierr<0) { + char * groupname = CCTK_GroupName(gi); + assert (groupname); + CCTK_VWarn (0, __LINE__, __FILE__, CCTK_THORNSTRING, + "Error in tensor type alias declaration for group \"%s\"", + groupname); + free (groupname); + } + + if (CCTK_EQUALS (tensortypealias, "scalar")) { + /* scalar */ + if (numvars != 1) { + char * groupname = CCTK_GroupName(gi); + assert (groupname); + CCTK_VWarn (2, __LINE__, __FILE__, CCTK_THORNSTRING, + "Group \"%s\" has the tensor type alias \"scalar\", but contains more than 1 element", + groupname); + free (groupname); + } + } else if (CCTK_EQUALS (tensortypealias, "u")) { + /* vector */ + assert (numvars == 3); + } else if (CCTK_EQUALS (tensortypealias, "dd_sym")) { + /* symmetric tensor */ + assert (numvars == 6); + } else { + char * groupname = CCTK_GroupName(gi); + assert (groupname); + CCTK_VWarn (0, __LINE__, __FILE__, CCTK_THORNSTRING, + "Illegal tensor type alias for group \"%s\"", + groupname); + free (groupname); + } + } + + { + int min_handle, max_handle; + CCTK_REAL local[6], global[6]; + min_handle = CCTK_ReductionArrayHandle ("minimum"); + if (min_handle<0) CCTK_WARN (0, "Could not obtain reduction handle"); + max_handle = CCTK_ReductionArrayHandle ("maximum"); + if (max_handle<0) CCTK_WARN (0, "Could not obtain reduction handle"); + + for (d=0; d<6; ++d) local[d] = cctkGH->cctk_bbox[d]; + ierr = CCTK_ReduceLocArrayToArray1D + (cctkGH, -1, max_handle, local, global, 6, CCTK_VARIABLE_REAL); + for (d=0; d<6; ++d) global_bbox[d] = (int)global[d]; + + for (d=0; d<3; ++d) local[d] = cctkGH->cctk_lbnd[d]; + ierr = CCTK_ReduceLocArrayToArray1D + (cctkGH, -1, min_handle, local, global, 3, CCTK_VARIABLE_REAL); + for (d=0; d<3; ++d) global_lbnd[d] = (int)global[d]; + + for (d=0; d<3; ++d) local[d] = cctkGH->cctk_ubnd[d]; + ierr = CCTK_ReduceLocArrayToArray1D + (cctkGH, -1, max_handle, local, global, 3, CCTK_VARIABLE_REAL); + for (d=0; d<3; ++d) global_ubnd[d] = (int)global[d]; + + for (d=0; d<3; ++d) { + fake_bbox[2*d ] = data.lbnd[d] == global_lbnd[d]; + fake_bbox[2*d+1] = data.ubnd[d] == global_ubnd[d]; + } + } + + /* directions */ + alldirs[0] = 0; + alldirs[1] = 1; + + /* Find grid point that corresponds to the origin */ + for (q=0; q<2; ++q) { + d = alldirs[q]; + /* x0 + dx * origin == 0 */ + origin[d] = - x0[d] / dx[d]; + dorigin[d] = origin[d] - floor(origin[d]); + if (fabs(dorigin[d]) < 1.0e-6 || fabs(dorigin[d] - 1.0) < 1.0e-6) { + avoid_origin[d] = 0; + } else if (fabs(dorigin[d] - 0.5) < 1.0e-6) { + avoid_origin[d] = 1; + } else { + CCTK_VWarn (0, __LINE__, __FILE__, CCTK_THORNSTRING, + "The coordinate origin in the %c-direction falls neither onto a grid point nor into the middle between two grid points.", + "xyz"[d]); + } + offset[d] = avoid_origin[d] ? 0 : 1; + } + + for (step=0; step<3; ++step) { + switch (step) { + case 0: + /* x face */ + ndirs = 1; + dir[0] = 0; + otherdir[0] = 1; + break; + case 1: + /* y face */ + ndirs = 1; + dir[0] = 1; + otherdir[0] = 0; + break; + case 2: + /* xy edge */ + ndirs = 2; + dir[0] = 0; + dir[1] = 1; + otherdir[0] = 1; + otherdir[1] = 0; + break; + default: + assert (0); + } + + if (CCTK_EQUALS (tensortypealias, "scalar")) { + /* scalar */ + srcvi = vi; + parity = +1; + } else if (CCTK_EQUALS (tensortypealias, "u")) { + /* vector */ + int srcindex; + if (srcindex == dir[0]) srcindex = otherdir[0], parity = +1; + else if (srcindex == otherdir[0]) srcindex = dir[0], parity = -1; + else srcindex = index, parity = +1; + srcvi = firstvar + srcindex; + } else if (CCTK_EQUALS (tensortypealias, "dd_sym")) { + /* symmetric tensor */ + int index1, index2; + int srcindex1, srcindex2; + int srcindex; + switch (index) { + case 0: index1=0, index2=0; break; + case 1: index1=0, index2=1; break; + case 2: index1=0, index2=2; break; + case 3: index1=1, index2=1; break; + case 4: index1=1, index2=2; break; + case 5: index1=2, index2=2; break; + default: assert(0); + } + parity = +1; + if (srcindex1 == dir[0]) srcindex1 = otherdir[0], parity *= +1; + else if (srcindex1 == otherdir[0]) srcindex1 = dir[0], parity *= -1; + else srcindex1 = index1, parity *= +1; + if (srcindex2 == dir[0]) srcindex2 = otherdir[0], parity *= +1; + else if (srcindex2 == otherdir[0]) srcindex2 = dir[0], parity *= -1; + else srcindex2 = index2, parity *= +1; + { + int compact[3][3] = { 0,1,2, 1,3,4, 2,4,5 }; + srcindex = compact[srcindex1][srcindex2]; + } + srcvi = firstvar + srcindex; + } else { + assert (0); + } + srcptr = CCTK_VarDataPtrI (cctkGH, 0, srcvi); + assert (srcptr); + + have_global_bbox = 1; + for (q=0; q= 0); + + if (data.gsh[otherdir[q]] < stencil[dir[q]] + stencil[otherdir[q]] + offset[otherdir[q]]) { + CCTK_VWarn (0, __LINE__, __FILE__, CCTK_THORNSTRING, + "The group \"%s\" has in the %c-direction only %d grid points. " + "This is not large enough for a 90 degree rotating symmetry boundary that is %d grid points wide in the %c-direction. " + "The group needs to have at least %d grid points in the %c-direction.", + CCTK_GroupNameFromVarI(vi), "xyz"[dir[q]], data.gsh[dir[q]], + stencil[dir[q]], "xyz"[otherdir[q]], + stencil[dir[q]] + stencil[otherdir[q]] + offset[otherdir[q]], "xyz"[dir[q]]); + } + } + + /* Allocate slab transfer description */ + xferinfo = malloc(group.dim * sizeof *xferinfo); + assert (xferinfo); + + /* Fill in slab transfer description */ + for (d=0; d=0); + + indices = malloc (nvars * sizeof *indices); + assert (indices); + faces = malloc (nvars * sizeof *faces); + assert (faces); + widths = malloc (nvars * sizeof *widths); + assert (widths); + tables = malloc (nvars * sizeof *tables); + assert (tables); + + ierr = Boundary_SelectedGVs + (cctkGH, nvars, indices, faces, widths, tables, 0); + assert (ierr == nvars); + + for (i=0; i=0 && vi= 0); + + dim = CCTK_GroupDimFromVarI (vi); + assert (dim>=0); + + stencil = malloc (dim * sizeof *stencil); + assert (stencil); + ierr = CCTK_GroupnghostzonesVI (cctkGH, dim, stencil, vi); + assert (!ierr); + + ierr = BndRot90VI (cctkGH, stencil, vi); + assert (!ierr); + + free (stencil); + } + + free (indices); + free (faces); + free (widths); + free (tables); +} -- cgit v1.2.3