Fixed Neutron Source Diffusion

18 views
Skip to first unread message

Joe Gayo

unread,
Jul 9, 2021, 3:02:52 PM7/9/21
to wasora
Hello,

As a starting point, I'm attempting to model neutron diffusion through a cube of water from a fixed source. The setup is a thin volume on one face of the main water cube acting as the source. It doesn't seem to model the diffusion between the two materials. I'm uncertain what boundary condition or setting needs to be implemented to enable this. I've attached an image of the result (the group doesn't seem to allow the posting of .mil or .geo files).

Thanks,

Joe
Boundary.PNG

Joe Gayo

unread,
Jul 9, 2021, 3:04:01 PM7/9/21
to wasora, Joe Gayo
# if the mesh does not already exists, call gmsh on the $1.geo file
SHELL "if [ ! -e cubeWater.msh ]; then gmsh -v 0 -3 cubeWater.geo; fi"
MESH FILE_PATH cubeWater.msh DIMENSIONS 3

# these are the default values, they are overwritten by the commandline arguments
MILONGA_PROBLEM GROUPS 1 SCHEME elements FORMULATION diffusion

# some settings to improve cpu & memory usage
MILONGA_SOLVER EPS_TYPE jd ST_TYPE precond KSP_TYPE bcgs PC_TYPE asm

# table 12 in page 18 of Los Alamos Report LA-13511
# analytical benchmark test set for criticality code verification
MATERIAL Water S 0 SigmaT 0.32640 SigmaS 0.293760 SigmaA 0.032640 nuSigmaF 0.0
MATERIAL Air   S 1 D 1

# link physical entities in the mesh to materials and boundary conditions
MATERIAL Water PHYSICAL_ENTITY moderator
MATERIAL Air PHYSICAL_ENTITY source
PHYSICAL_ENTITY external BC vacuum
PHYSICAL_ENTITY source_mirror BC mirror
#PHYSICAL_ENTITY interface


# sn_alpha = 1.0

# do the magic!
MILONGA_STEP

# write some results into to the standard output
# PRINT_FUNCTION phi1 HEADER

Joe Gayo

unread,
Jul 9, 2021, 3:05:18 PM7/9/21
to wasora, Joe Gayo
//
SetFactory("OpenCASCADE");


a = 12.5;        // cube half-length (cm)
b = 1.0;          // source sheet
lc = a/12.5;        // element characteristic length
Mesh.CharacteristicLengthMin = lc;
Mesh.CharacteristicLengthMax = lc;

Box (1) = {-(a+b), -a, -a, b, 2*a, 2*a};
Box (2) = {-a, -a, -a, 2*a, 2*a, 2*a};

Physical Surface("external") = {3,4,5,6,9,10,11,12,8};
Physical Surface("source_mirror") = {1};
Physical Surface("interface") = {2, 7};
Physical Volume("source") = {1};
Physical Volume("moderator") = {2};


Mesh.Optimize = 1;
Mesh.OptimizeNetgen = 0;

Mesh.Algorithm = 6;    
// (1=MeshAdapt, 2=Automatic, 5=Delaunay, 6=Frontal, 7=BAMG, 8=DelQuad)
Mesh.Algorithm3D = 4;  
// (1=Delaunay, 4=Frontal, 5=Frontal Delaunay, 6=Frontal Hex, 7=MMG3D, 9=R-tree)

Joe Gayo

unread,
Jul 9, 2021, 3:52:03 PM7/9/21
to wasora
I used the "Coherence;" command in the .geo file and renumbered the surfaces. This seemed to work. Maybe there was a double surface at the interface causing the problem. I still don't know if I'm doing this correctly (the new plot could be correct). 

Jeremy Theler

unread,
Jul 9, 2021, 6:31:17 PM7/9/21
to was...@seamplex.com
Indeed, when you create two boxes eache of them has a unique set of 6 boundary faces.
When you call Coherence (or Boolean Fragments) then the common faces are identified and joined together.

Ironically, I cannot compile milonga right now so  I cannot test your case, but

 1. you don't need to set a physical group for the interface
 2. if you used two groups the physics would be more interesting
 3. if you set the "external" surface to have a mirror condition you can compare the results with a 1d slab (which I believe has analytical solution)

You can also use the attached geo file that creates a structured hexahedral mesh.

Let us know how it goes
--
jeremy theler
--
You received this message because you are subscribed to the Google Groups "wasora" group.
To unsubscribe from this group and stop receiving emails from it, send an email to wasora+un...@seamplex.com.
To view this discussion on the web visit https://groups.google.com/a/seamplex.com/d/msgid/wasora/1bd02c85-8361-4966-b1a4-7333f4d7e439n%40seamplex.com.

joe-struct.msh

Jeremy Theler

unread,
Jul 9, 2021, 7:02:59 PM7/9/21
to was...@seamplex.com
I had milonga run on a debian buster VM becasue it no longer compiles with gcc >= 9. I am now working on a new tool that will include wasora, milonga and Fino at

Here's the result I got with the structured mesh. It looks good, doesn't it?


On Fri, 2021-07-09 at 12:52 -0700, Joe Gayo wrote:
--
joe.png

Joe Gayo

unread,
Jul 9, 2021, 7:09:03 PM7/9/21
to wasora, Jeremy
Thanks Jeremy.

I haven't run the structured mesh yet, but your image looks good.

Attached is the updated image from my original mesh.

I'm having a hard time with the 1D solution:

MESH FILE_PATH reedWater.msh
MILONGA_PROBLEM DIMENSIONS 1 GROUPS 1 FORMULATION diffusion
sn_alpha = 1.0

MATERIAL source        S 1  D 1
#Table 12 in page 18 of Los Alamos Report LA-13511, Water
MATERIAL moderator     S 0 SigmaT 0.32640 SigmaS 0.293760 SigmaA 0.032640 nuSigmaF 0.0

#If there physical entity named "external" or "mirror" gets the following BC
PHYSICAL_ENTITY external BC vacuum
PHYSICAL_ENTITY source_mirror BC mirror

MILONGA_STEP

PRINT_FUNCTION phi1 HEADER



lc = 0.01;

Point(1) = {0, 0, 0, lc};
Point(2) = {0.1, 0, 0, lc};
Point(3) = {12.6, 0, 0, lc};

Point(11) = {0.05, 0, 0, lc};
Point(12) = {6.35, 0, 0, lc};

Line(1) = {1, 11};
Line(11) = {11, 2};
Line(2) = {2, 12};
Line(12) = {12, 3};

Physical Line("source")  = {1,11};
Physical Line("moderator") = {2,12};

Physical Point("source_mirror") = {1};
Physical Point("external") = {3};

WaterCube.PNG

Joe Gayo

unread,
Jul 9, 2021, 7:22:28 PM7/9/21
to wasora, Joe Gayo, Jeremy
Added this line from an example file...

# some settings to improve cpu & memory usage
MILONGA_SOLVER EPS_TYPE jd ST_TYPE precond KSP_TYPE bcgs PC_TYPE asm

and the 1D now solves.

The structured mesh is great! How do I create that?

Jeremy Theler

unread,
Jul 9, 2021, 7:53:06 PM7/9/21
to Joe Gayo, wasora
What was the error you were getting for the 1d case?

It is better if you attach files instead of pasting the contents so we can track changes.
I attached the geo file to create the structured mesh in a previous email, it was called joe-struct.geo
The last png you shared looks great. Do you have XSs for two groups?

Joe Gayo

unread,
Jul 9, 2021, 7:59:54 PM7/9/21
to wasora, Jeremy, Joe Gayo
I don't see the geo file.

I was pasting the code because it wasn't allowing me to post with geo and mil files attached.

Here is the error:
gayo8875@gayo8875-VirtualBox:~/milonga/examples$ ./test-reedWater.sh
error: PETSc's linear solver did not converge with reason 'DIVERGED_ITS' (-3)
"reedWater.gp" line 6: warning: Skipping data file with no valid points

plot 'reedWater.dat'   w p pt 6 lt 1,
                                     ^
"reedWater.gp" line 6: all points y value undefined!


I need to find XSs for water and paraffin (or Polyethylene) for two groups.

Jeremy Theler

unread,
Jul 9, 2021, 8:06:09 PM7/9/21
to Joe Gayo, wasora
Here's the geo file again.
Weird that it does not converge, but this is what I said before that I did not tweak the linear solvers for the fixed-source case. The solvers are optimized for the multiplicative media without sources (i.e. computation of keff).
It is also weird that you were not allowed to add attachments.

Regards
--
jeremy
joe-struct.geo

Joe Gayo

unread,
Jul 9, 2021, 8:08:34 PM7/9/21
to wasora, Jeremy, Joe Gayo
Do you know where I can finds XSs for water and paraffin (or Polyethylene) for two groups?

Joe Gayo

unread,
Jul 9, 2021, 10:16:26 PM7/9/21
to wasora, Joe Gayo, Jeremy
Is it possible to obtain a plot of fluence for neutrons within a certain energy band? For example, fast neutrons enter and where is the high fluence of thermal neutrons in the moderator?
Reply all
Reply to author
Forward
0 new messages