1 static char help[] = "Vlasov example of particles orbiting around a central massive point.\n"; 2 3 #include <petscdmplex.h> 4 #include <petsc/private/dmpleximpl.h> /* For norm */ 5 #include <petsc/private/petscfeimpl.h> /* For CoordinatesRefToReal() */ 6 #include <petscdmswarm.h> 7 #include <petscts.h> 8 9 typedef struct { 10 PetscInt dim; /* The topological mesh dimension */ 11 PetscBool simplex; /* Flag for simplices or tensor cells */ 12 char filename[PETSC_MAX_PATH_LEN]; /* Name of the mesh filename if any */ 13 PetscInt particlesPerCell; /* The number of partices per cell */ 14 PetscReal momentTol; /* Tolerance for checking moment conservation */ 15 PetscBool monitor; /* Flag for using the TS monitor */ 16 PetscBool error; /* Flag for printing the error */ 17 PetscInt ostep; /* print the energy at each ostep time steps */ 18 } AppCtx; 19 20 static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 21 { 22 PetscErrorCode ierr; 23 24 PetscFunctionBeginUser; 25 options->dim = 2; 26 options->simplex = PETSC_TRUE; 27 options->monitor = PETSC_FALSE; 28 options->error = PETSC_FALSE; 29 options->particlesPerCell = 1; 30 options->momentTol = 100.0*PETSC_MACHINE_EPSILON; 31 options->ostep = 100; 32 33 ierr = PetscStrcpy(options->filename, "");CHKERRQ(ierr); 34 35 ierr = PetscOptionsBegin(comm, "", "Vlasov Options", "DMPLEX");CHKERRQ(ierr); 36 ierr = PetscOptionsInt("-output_step", "Number of time steps between output", "ex4.c", options->ostep, &options->ostep, PETSC_NULL);CHKERRQ(ierr); 37 ierr = PetscOptionsInt("-dim", "The topological mesh dimension", "ex4.c", options->dim, &options->dim, NULL);CHKERRQ(ierr); 38 ierr = PetscOptionsBool("-monitor", "Flag to use the TS monitor", "ex4.c", options->monitor, &options->monitor, NULL);CHKERRQ(ierr); 39 ierr = PetscOptionsBool("-error", "Flag to print the error", "ex4.c", options->error, &options->error, NULL);CHKERRQ(ierr); 40 ierr = PetscOptionsBool("-simplex", "The flag for simplices or tensor cells", "ex4.c", options->simplex, &options->simplex, NULL);CHKERRQ(ierr); 41 ierr = PetscOptionsString("-mesh", "Name of the mesh filename if any", "ex4.c", options->filename, options->filename, PETSC_MAX_PATH_LEN, NULL);CHKERRQ(ierr); 42 ierr = PetscOptionsInt("-particles_per_cell", "Number of particles per cell", "ex4.c", options->particlesPerCell, &options->particlesPerCell, NULL);CHKERRQ(ierr); 43 ierr = PetscOptionsEnd();CHKERRQ(ierr); 44 45 PetscFunctionReturn(0); 46 } 47 48 static PetscErrorCode CreateMesh(MPI_Comm comm, DM *dm, AppCtx *user) 49 { 50 PetscBool flg; 51 PetscErrorCode ierr; 52 53 PetscFunctionBeginUser; 54 ierr = PetscStrcmp(user->filename, "", &flg);CHKERRQ(ierr); 55 if (flg) { 56 ierr = DMPlexCreateBoxMesh(comm, user->dim, user->simplex, NULL, NULL, NULL, NULL, PETSC_TRUE, dm);CHKERRQ(ierr); 57 } else { 58 ierr = DMPlexCreateFromFile(comm, user->filename, PETSC_TRUE, dm);CHKERRQ(ierr); 59 ierr = DMGetDimension(*dm, &user->dim);CHKERRQ(ierr); 60 } 61 { 62 DM distributedMesh = NULL; 63 64 ierr = DMPlexDistribute(*dm, 0, NULL, &distributedMesh);CHKERRQ(ierr); 65 if (distributedMesh) { 66 ierr = DMDestroy(dm);CHKERRQ(ierr); 67 *dm = distributedMesh; 68 } 69 } 70 ierr = DMLocalizeCoordinates(*dm);CHKERRQ(ierr); /* needed for periodic */ 71 ierr = DMSetFromOptions(*dm);CHKERRQ(ierr); 72 ierr = PetscObjectSetName((PetscObject) *dm, "Mesh");CHKERRQ(ierr); 73 ierr = DMViewFromOptions(*dm, NULL, "-dm_view");CHKERRQ(ierr); 74 PetscFunctionReturn(0); 75 } 76 77 static PetscErrorCode SetInitialCoordinates(DM dmSw) 78 { 79 DM dm; 80 AppCtx *user; 81 PetscRandom rnd; 82 PetscBool simplex; 83 PetscReal *centroid, *coords, *xi0, *v0, *J, *invJ, detJ; 84 PetscInt dim, d, cStart, cEnd, c, Np, p; 85 PetscErrorCode ierr; 86 87 PetscFunctionBeginUser; 88 ierr = PetscRandomCreate(PetscObjectComm((PetscObject) dmSw), &rnd);CHKERRQ(ierr); 89 ierr = PetscRandomSetInterval(rnd, -1.0, 1.0);CHKERRQ(ierr); 90 ierr = PetscRandomSetFromOptions(rnd);CHKERRQ(ierr); 91 92 ierr = DMGetApplicationContext(dmSw, (void **) &user);CHKERRQ(ierr); 93 simplex = user->simplex; 94 Np = user->particlesPerCell; 95 ierr = DMSwarmGetCellDM(dmSw, &dm);CHKERRQ(ierr); 96 ierr = DMGetDimension(dm, &dim);CHKERRQ(ierr); 97 ierr = DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd);CHKERRQ(ierr); 98 ierr = PetscMalloc5(dim, ¢roid, dim, &xi0, dim, &v0, dim*dim, &J, dim*dim, &invJ);CHKERRQ(ierr); 99 for (d = 0; d < dim; ++d) xi0[d] = -1.0; 100 ierr = DMSwarmGetField(dmSw, DMSwarmPICField_coor, NULL, NULL, (void **) &coords);CHKERRQ(ierr); 101 for (c = cStart; c < cEnd; ++c) { 102 if (Np == 1) { 103 ierr = DMPlexComputeCellGeometryFVM(dm, c, NULL, centroid, NULL);CHKERRQ(ierr); 104 for (d = 0; d < dim; ++d) coords[c*dim+d] = centroid[d]; 105 } else { 106 ierr = DMPlexComputeCellGeometryFEM(dm, c, NULL, v0, J, invJ, &detJ);CHKERRQ(ierr); /* affine */ 107 for (p = 0; p < Np; ++p) { 108 const PetscInt n = c*Np + p; 109 PetscReal sum = 0.0, refcoords[3]; 110 111 for (d = 0; d < dim; ++d) { 112 ierr = PetscRandomGetValueReal(rnd, &refcoords[d]);CHKERRQ(ierr); 113 sum += refcoords[d]; 114 } 115 if (simplex && sum > 0.0) for (d = 0; d < dim; ++d) refcoords[d] -= PetscSqrtReal(dim)*sum; 116 CoordinatesRefToReal(dim, dim, xi0, v0, J, refcoords, &coords[n*dim]); 117 } 118 } 119 } 120 ierr = DMSwarmRestoreField(dmSw, DMSwarmPICField_coor, NULL, NULL, (void **) &coords);CHKERRQ(ierr); 121 ierr = PetscFree5(centroid, xi0, v0, J, invJ);CHKERRQ(ierr); 122 ierr = PetscRandomDestroy(&rnd);CHKERRQ(ierr); 123 PetscFunctionReturn(0); 124 } 125 126 static PetscErrorCode SetInitialConditions(DM dmSw, Vec u) 127 { 128 DM dm; 129 AppCtx *user; 130 PetscScalar *initialConditions; 131 PetscInt dim, cStart, cEnd, c, Np, p; 132 PetscErrorCode ierr; 133 134 PetscFunctionBeginUser; 135 ierr = DMGetApplicationContext(dmSw, (void **) &user);CHKERRQ(ierr); 136 Np = user->particlesPerCell; 137 ierr = DMSwarmGetCellDM(dmSw, &dm);CHKERRQ(ierr); 138 ierr = DMGetDimension(dm, &dim);CHKERRQ(ierr); 139 ierr = DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd);CHKERRQ(ierr); 140 ierr = VecGetArray(u, &initialConditions);CHKERRQ(ierr); 141 for (c = cStart; c < cEnd; ++c) { 142 for (p = 0; p < Np; ++p) { 143 const PetscInt n = c*Np + p; 144 145 initialConditions[(n*2 + 0)*dim + 0] = n+1; 146 initialConditions[(n*2 + 0)*dim + 1] = 0; 147 initialConditions[(n*2 + 1)*dim + 0] = 0; 148 initialConditions[(n*2 + 1)*dim + 1] = PetscSqrtReal(1000./(n+1.)); 149 } 150 } 151 ierr = VecRestoreArray(u, &initialConditions);CHKERRQ(ierr); 152 PetscFunctionReturn(0); 153 } 154 155 static PetscErrorCode CreateParticles(DM dm, DM *sw, AppCtx *user) 156 { 157 PetscInt *cellid; 158 PetscInt dim, cStart, cEnd, c, Np = user->particlesPerCell, p; 159 PetscErrorCode ierr; 160 161 PetscFunctionBeginUser; 162 ierr = DMGetDimension(dm, &dim);CHKERRQ(ierr); 163 ierr = DMCreate(PetscObjectComm((PetscObject) dm), sw);CHKERRQ(ierr); 164 ierr = DMSetType(*sw, DMSWARM);CHKERRQ(ierr); 165 ierr = DMSetDimension(*sw, dim);CHKERRQ(ierr); 166 167 ierr = DMSwarmSetType(*sw, DMSWARM_PIC);CHKERRQ(ierr); 168 ierr = DMSwarmSetCellDM(*sw, dm);CHKERRQ(ierr); 169 ierr = DMSwarmRegisterPetscDatatypeField(*sw, "kinematics", 2*dim, PETSC_REAL);CHKERRQ(ierr); 170 ierr = DMSwarmFinalizeFieldRegister(*sw);CHKERRQ(ierr); 171 ierr = DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd);CHKERRQ(ierr); 172 ierr = DMSwarmSetLocalSizes(*sw, (cEnd - cStart) * Np, 0);CHKERRQ(ierr); 173 ierr = DMSetFromOptions(*sw);CHKERRQ(ierr); 174 ierr = DMSwarmGetField(*sw, DMSwarmPICField_cellid, NULL, NULL, (void **) &cellid);CHKERRQ(ierr); 175 for (c = cStart; c < cEnd; ++c) { 176 for (p = 0; p < Np; ++p) { 177 const PetscInt n = c*Np + p; 178 179 cellid[n] = c; 180 } 181 } 182 ierr = DMSwarmRestoreField(*sw, DMSwarmPICField_cellid, NULL, NULL, (void **) &cellid);CHKERRQ(ierr); 183 ierr = PetscObjectSetName((PetscObject) *sw, "Particles");CHKERRQ(ierr); 184 ierr = DMViewFromOptions(*sw, NULL, "-sw_view");CHKERRQ(ierr); 185 PetscFunctionReturn(0); 186 } 187 188 /* Create particle RHS Functions for gravity with G = 1 for simplification */ 189 static PetscErrorCode RHSFunction1(TS ts, PetscReal t, Vec V, Vec Xres, void *ctx) 190 { 191 DM dm; 192 const PetscScalar *v; 193 PetscScalar *xres; 194 PetscInt Np, p, dim, d; 195 PetscErrorCode ierr; 196 197 PetscFunctionBeginUser; 198 ierr = TSGetDM(ts, &dm);CHKERRQ(ierr); 199 ierr = DMGetDimension(dm, &dim);CHKERRQ(ierr); 200 ierr = VecGetLocalSize(Xres, &Np);CHKERRQ(ierr); 201 Np /= dim; 202 ierr = VecGetArray(Xres, &xres);CHKERRQ(ierr); 203 ierr = VecGetArrayRead(V, &v);CHKERRQ(ierr); 204 for (p = 0; p < Np; ++p) { 205 for(d = 0; d < dim; ++d) { 206 xres[p*dim+d] = v[p*dim+d]; 207 } 208 } 209 ierr = VecRestoreArrayRead(V,& v);CHKERRQ(ierr); 210 ierr = VecRestoreArray(Xres, &xres);CHKERRQ(ierr); 211 PetscFunctionReturn(0); 212 } 213 214 static PetscErrorCode RHSFunction2(TS ts, PetscReal t, Vec X, Vec Vres, void *ctx) 215 { 216 DM dm; 217 const PetscScalar *x; 218 PetscScalar *vres; 219 PetscInt Np, p, dim, d; 220 PetscErrorCode ierr; 221 222 223 PetscFunctionBeginUser; 224 ierr = TSGetDM(ts, &dm);CHKERRQ(ierr); 225 ierr = DMGetDimension(dm, &dim);CHKERRQ(ierr); 226 ierr = VecGetLocalSize(Vres, &Np);CHKERRQ(ierr); 227 Np /= dim; 228 ierr = VecGetArray(Vres, &vres);CHKERRQ(ierr); 229 ierr = VecGetArrayRead(X, &x);CHKERRQ(ierr); 230 for (p = 0; p < Np; ++p) { 231 const PetscScalar rsqr = DMPlex_NormD_Internal(dim, &x[p*dim]); 232 233 for (d = 0; d < dim; ++d) { 234 vres[p*dim+d] = -(1000./(p+1.)) * x[p*dim+d]/PetscSqr(rsqr); 235 } 236 } 237 ierr = VecRestoreArray(Vres, &vres);CHKERRQ(ierr); 238 ierr = VecRestoreArrayRead(X, &x);CHKERRQ(ierr); 239 PetscFunctionReturn(0); 240 } 241 242 static PetscErrorCode RHSFunctionParticles(TS ts, PetscReal t , Vec U, Vec R, void *ctx) 243 { 244 DM dm; 245 const PetscScalar *u; 246 PetscScalar *r; 247 PetscInt Np, p, dim, d; 248 PetscErrorCode ierr; 249 250 PetscFunctionBeginUser; 251 ierr = TSGetDM(ts, &dm);CHKERRQ(ierr); 252 ierr = DMGetDimension(dm, &dim);CHKERRQ(ierr); 253 ierr = VecGetLocalSize(U, &Np);CHKERRQ(ierr); 254 Np /= 2*dim; 255 ierr = VecGetArrayRead(U, &u);CHKERRQ(ierr); 256 ierr = VecGetArray(R, &r);CHKERRQ(ierr); 257 for( p = 0; p < Np; ++p) { 258 const PetscScalar rsqr = DMPlex_NormD_Internal(dim, &u[p*2*dim]); 259 260 for (d = 0; d < dim; ++d) { 261 r[p*2*dim+d] = u[p*2*dim+d+2]; 262 r[p*2*dim+d+2] = -(1000./(1.+p)) * u[p*2*dim+d]/PetscSqr(rsqr); 263 } 264 } 265 ierr = VecRestoreArrayRead(U, &u);CHKERRQ(ierr); 266 ierr = VecRestoreArray(R, &r);CHKERRQ(ierr); 267 PetscFunctionReturn(0); 268 } 269 270 static PetscErrorCode InitializeSolve(TS ts, Vec u) 271 { 272 DM dm; 273 AppCtx *user; 274 PetscErrorCode ierr; 275 276 PetscFunctionBeginUser; 277 ierr = TSGetDM(ts, &dm);CHKERRQ(ierr); 278 ierr = DMGetApplicationContext(dm, (void **) &user);CHKERRQ(ierr); 279 ierr = SetInitialCoordinates(dm);CHKERRQ(ierr); 280 ierr = SetInitialConditions(dm, u);CHKERRQ(ierr); 281 PetscFunctionReturn(0); 282 } 283 284 static PetscErrorCode ComputeError(TS ts, Vec U, Vec E) 285 { 286 MPI_Comm comm; 287 DM sdm; 288 AppCtx *user; 289 const PetscScalar *u, *coords; 290 PetscScalar *e; 291 PetscReal t; 292 PetscInt dim, Np, p; 293 PetscErrorCode ierr; 294 295 PetscFunctionBeginUser; 296 ierr = PetscObjectGetComm((PetscObject) ts, &comm);CHKERRQ(ierr); 297 ierr = TSGetDM(ts, &sdm);CHKERRQ(ierr); 298 ierr = DMGetApplicationContext(sdm, (void **) &user);CHKERRQ(ierr); 299 ierr = DMGetDimension(sdm, &dim);CHKERRQ(ierr); 300 ierr = TSGetSolveTime(ts, &t);CHKERRQ(ierr); 301 ierr = VecGetArray(E, &e);CHKERRQ(ierr); 302 ierr = VecGetArrayRead(U, &u);CHKERRQ(ierr); 303 ierr = VecGetLocalSize(U, &Np);CHKERRQ(ierr); 304 Np /= 2*dim; 305 ierr = DMSwarmGetField(sdm, DMSwarmPICField_coor, NULL, NULL, (void **) &coords);CHKERRQ(ierr); 306 for (p = 0; p < Np; ++p) { 307 const PetscScalar *x = &u[(p*2+0)*dim]; 308 const PetscScalar *v = &u[(p*2+1)*dim]; 309 const PetscReal x0 = p+1.; 310 const PetscReal omega = PetscSqrtReal(1000./(p+1.))/x0; 311 const PetscReal xe[3] = { x0*PetscCosReal(omega*t), x0*PetscSinReal(omega*t), 0.0}; 312 const PetscReal ve[3] = {-x0*omega*PetscSinReal(omega*t), x0*omega*PetscCosReal(omega*t), 0.0}; 313 PetscInt d; 314 315 for (d = 0; d < dim; ++d) { 316 e[(p*2+0)*dim+d] = x[d] - xe[d]; 317 e[(p*2+1)*dim+d] = v[d] - ve[d]; 318 } 319 if (user->error) {ierr = PetscPrintf(comm, "p%D error [%.2g %.2g] sol [(%.6lf %.6lf) (%.6lf %.6lf)] exact [(%.6lf %.6lf) (%.6lf %.6lf)] energy/exact energy %g / %g\n", p, (double) DMPlex_NormD_Internal(dim, &e[(p*2+0)*dim]), (double) DMPlex_NormD_Internal(dim, &e[(p*2+1)*dim]), (double) x[0], (double) x[1], (double) v[0], (double) v[1], (double) xe[0], (double) xe[1], (double) ve[0], (double) ve[1], 0.5*DMPlex_NormD_Internal(dim, v), (double) (0.5*(1000./(p+1))));} 320 } 321 ierr = DMSwarmRestoreField(sdm, DMSwarmPICField_coor, NULL, NULL, (void **) &coords);CHKERRQ(ierr); 322 ierr = VecRestoreArrayRead(U, &u);CHKERRQ(ierr); 323 ierr = VecRestoreArray(E, &e);CHKERRQ(ierr); 324 PetscFunctionReturn(0); 325 } 326 327 int main(int argc, char **argv) 328 { 329 TS ts; 330 TSConvergedReason reason; 331 DM dm, sw; 332 Vec u; 333 IS is1, is2; 334 PetscInt *idx1, *idx2; 335 MPI_Comm comm; 336 AppCtx user; 337 const PetscScalar *endVals; 338 PetscReal ftime = .1; 339 PetscInt locSize, dim, d, Np, p, steps; 340 PetscErrorCode ierr; 341 342 ierr = PetscInitialize(&argc, &argv, NULL, help);if (ierr) return ierr; 343 comm = PETSC_COMM_WORLD; 344 ierr = ProcessOptions(comm, &user);CHKERRQ(ierr); 345 346 ierr = CreateMesh(comm, &dm, &user);CHKERRQ(ierr); 347 ierr = CreateParticles(dm, &sw, &user);CHKERRQ(ierr); 348 ierr = DMSetApplicationContext(sw, &user);CHKERRQ(ierr); 349 350 ierr = TSCreate(comm, &ts);CHKERRQ(ierr); 351 ierr = TSSetType(ts, TSBASICSYMPLECTIC);CHKERRQ(ierr); 352 ierr = TSSetDM(ts, sw);CHKERRQ(ierr); 353 ierr = TSSetMaxTime(ts, ftime);CHKERRQ(ierr); 354 ierr = TSSetTimeStep(ts, 0.0001);CHKERRQ(ierr); 355 ierr = TSSetMaxSteps(ts, 10);CHKERRQ(ierr); 356 ierr = TSSetExactFinalTime(ts, TS_EXACTFINALTIME_MATCHSTEP);CHKERRQ(ierr); 357 ierr = TSSetTime(ts, 0.0);CHKERRQ(ierr); 358 ierr = TSSetRHSFunction(ts, NULL, RHSFunctionParticles, &user);CHKERRQ(ierr); 359 360 ierr = DMSwarmCreateGlobalVectorFromField(sw, "kinematics", &u);CHKERRQ(ierr); 361 ierr = DMGetDimension(sw, &dim);CHKERRQ(ierr); 362 ierr = VecGetLocalSize(u, &locSize);CHKERRQ(ierr); 363 Np = locSize/(2*dim); 364 ierr = PetscMalloc1(locSize/2, &idx1);CHKERRQ(ierr); 365 ierr = PetscMalloc1(locSize/2, &idx2);CHKERRQ(ierr); 366 for (p = 0; p < Np; ++p) { 367 for (d = 0; d < dim; ++d) { 368 idx1[p*dim+d] = (p*2+0)*dim + d; 369 idx2[p*dim+d] = (p*2+1)*dim + d; 370 } 371 } 372 ierr = ISCreateGeneral(comm, locSize/2, idx1, PETSC_OWN_POINTER, &is1);CHKERRQ(ierr); 373 ierr = ISCreateGeneral(comm, locSize/2, idx2, PETSC_OWN_POINTER, &is2);CHKERRQ(ierr); 374 ierr = TSRHSSplitSetIS(ts, "position", is1);CHKERRQ(ierr); 375 ierr = TSRHSSplitSetIS(ts, "momentum", is2);CHKERRQ(ierr); 376 ierr = ISDestroy(&is1);CHKERRQ(ierr); 377 ierr = ISDestroy(&is2);CHKERRQ(ierr); 378 ierr = TSRHSSplitSetRHSFunction(ts,"position",NULL,RHSFunction1,&user);CHKERRQ(ierr); 379 ierr = TSRHSSplitSetRHSFunction(ts,"momentum",NULL,RHSFunction2,&user);CHKERRQ(ierr); 380 381 ierr = TSSetFromOptions(ts);CHKERRQ(ierr); 382 ierr = TSSetComputeInitialCondition(ts, InitializeSolve);CHKERRQ(ierr); 383 ierr = TSSetComputeExactError(ts, ComputeError);CHKERRQ(ierr); 384 ierr = TSComputeInitialCondition(ts, u);CHKERRQ(ierr); 385 ierr = VecViewFromOptions(u, NULL, "-init_view");CHKERRQ(ierr); 386 ierr = TSSolve(ts, u);CHKERRQ(ierr); 387 ierr = TSGetSolveTime(ts, &ftime);CHKERRQ(ierr); 388 ierr = TSGetConvergedReason(ts, &reason);CHKERRQ(ierr); 389 ierr = TSGetStepNumber(ts, &steps);CHKERRQ(ierr); 390 ierr = PetscPrintf(comm,"%s at time %g after %D steps\n", TSConvergedReasons[reason], (double) ftime, steps);CHKERRQ(ierr); 391 392 ierr = VecGetArrayRead(u, &endVals);CHKERRQ(ierr); 393 for (p = 0; p < Np; ++p) { 394 const PetscReal norm = DMPlex_NormD_Internal(dim, &endVals[(p*2 + 1)*dim]); 395 ierr = PetscPrintf(comm, "Particle %D initial Energy: %g Final Energy: %g\n", p, (double) (0.5*(1000./(p+1))), (double) 0.5*PetscSqr(norm));CHKERRQ(ierr); 396 } 397 ierr = VecRestoreArrayRead(u, &endVals);CHKERRQ(ierr); 398 ierr = DMSwarmDestroyGlobalVectorFromField(sw, "kinematics", &u);CHKERRQ(ierr); 399 ierr = TSDestroy(&ts);CHKERRQ(ierr); 400 ierr = DMDestroy(&sw);CHKERRQ(ierr); 401 ierr = DMDestroy(&dm);CHKERRQ(ierr); 402 ierr = PetscFinalize(); 403 return ierr; 404 } 405 406 407 /*TEST 408 409 build: 410 requires: triangle !single !complex 411 test: 412 suffix: bsi1 413 args: -dim 2 -dm_plex_box_faces 1,1 -dm_view -sw_view -ts_basicsymplectic_type 1 -ts_max_time 0.1 -ts_monitor_sp_swarm -ts_convergence_estimate -convest_num_refine 2 414 test: 415 suffix: bsi2 416 args: -dim 2 -dm_plex_box_faces 1,1 -dm_view -sw_view -ts_basicsymplectic_type 2 -ts_max_time 0.1 -ts_monitor_sp_swarm -ts_convergence_estimate -convest_num_refine 2 417 test: 418 suffix: euler 419 args: -dim 2 -dm_plex_box_faces 1,1 -dm_view -sw_view -ts_type euler -ts_max_time 0.1 -ts_monitor_sp_swarm -ts_convergence_estimate -convest_num_refine 2 420 421 TEST*/ 422