Vlasov Solver [work in progress]

Discuss how polywell fusion works; share theoretical questions and answers.

Moderators: tonybarry, MSimon

kcdodd
Posts: 722
Joined: Tue Jun 03, 2008 3:36 am
Location: Austin, TX

Post by kcdodd »

Well, if someone wanted to take a look I would not mind. But I know I hate reading other people's code trying to figure out what it's supposed to do, much less figure out why it's wrong. And it's not a small project. And it's in C++.
Carter

drmike
Posts: 825
Joined: Sat Jul 14, 2007 11:54 pm
Contact:

Post by drmike »

It's a hard problem, no matter how you do it! Take smaller bites to chew on and make things work one step at a time.

For my code, I'm just trying to get one species to work. It still doesn't work right. It's a fun challenge. Try to isolate problems to the point where you know something simple, and make it work correctly. You'll be back to that chunk of code again, no chunk is bug free. But at least it will be close to right so you can go to the next step.

I like having the whole frame work in place before I start debugging, but it's kind of silly. I end up commenting out most everything, debug a chunk, then uncomment some other chunk and debug that. Takes longer, but at least I see where I'm going!

Good luck. And don't give up!

kcdodd
Posts: 722
Joined: Tue Jun 03, 2008 3:36 am
Location: Austin, TX

Post by kcdodd »

Well,my biggest problems I think are coming from the fact that I don't know 6 dimensional geometry and I can't even really "see" what it's doing. I'm a very visual person I've just tried to piece together what I think it should/is doing. I have thought about rewriting the program work in lower number of dimensions so I could see it work but I would have to change a lot of code now.
Carter

drmike
Posts: 825
Joined: Sat Jul 14, 2007 11:54 pm
Contact:

Post by drmike »

No, you can leave all the code. Just don't loop over space dimensions. Pick one place (or no more than 10) and loop over velocity. If the 10 places are along a line, you have a 1D space problem and 3D (or 2D is better) in velocity. All the code stays - you just force the loops to be small or single valued along one dimension.

Then you can slowly unfold more dimensions and fix things which appear at each new level of complexity. But all the code basically stays in place.

kcdodd
Posts: 722
Joined: Tue Jun 03, 2008 3:36 am
Location: Austin, TX

Post by kcdodd »

Good news everyone! I love futurerama. I found one problem, yay.

p.s. boost::pool is worthless for large amounts of data.

I am going to make the project open source so people can steal my code for the betterment of mankind, haha. No, really my program is not very conventional. It is hard to explain the data structure, but in simplest terms the whole space is a binary partition tree. The lowest levels of the tree's partitions forms 6-D simplices, and the intersections of the partitions create the vertices. To integrate the density, a set of 3-D simplices (ie tetrahedrons) is created that encompass all momentum values. Those are then subdivided along the lines of the binary partition tree until the lowest level is reached where the average density in that 6-D simplex is averaged and multiplied by the volume of the resulting 3-D simplex, and then they are summed up to give the integral of a single point. I do this because the mesh is not a grid! It is irregular. Hopefully it will be on sourceforge soon.

However, I am starting to think it might be better to just make it a regular grid! The bad of an adaptive mesh seem to be overpowering the good.
Carter

tonybarry
Posts: 219
Joined: Sun Jul 08, 2007 4:32 am
Location: Sydney, Australia
Contact:

Post by tonybarry »

Hello Carter,
Thank you for deciding to make your code open source. Hopefully it will allow better bug detection.

If you can include as many details of your development environment that would be very helpful.
- Libraries (location and contents)
- makefile
- how the output is organised

etc etc

Regards,
Tony Barry

kcdodd
Posts: 722
Joined: Tue Jun 03, 2008 3:36 am
Location: Austin, TX

Post by kcdodd »

I created a project here: http://sourceforge.net/projects/avmi

I'm new to this stuff, especially cvs, so its a bit rough. The project is divided into three parts. DiscreteSpace is just a set of templates I made which allows the creation of an N dimensional space with N-simplex (triangle for N=2, tetrahedron for N=3, ...) subdivisions. VlasovSolver is the library that does all the simulations. TestApp is just a small console program that instantiates it all to test the functionality. The Util folder has some misc multiuse classes. I will have to work more on the documentation but the code is commented for the most part I think.

The dependencies are the 'matrix template library' and 'boost'. I don't have any makefiles, just codebocks project files.

And there is a bunch of garbage files (like old output and profiling data) I uploaded by accident. I still haven't figured out how to delete stuff.
Carter

kcdodd
Posts: 722
Joined: Tue Jun 03, 2008 3:36 am
Location: Austin, TX

Post by kcdodd »

When I run a test simulation it seems to exhibit undefined behaviour. One possibility I am thinking is that it is actually working but the integration routine is still buggy so it looks like it's not working. The other possibility is the stepping routine itself is unstable. I have an idea of a new stepping algorithm that might be more stable but if that is not even the problem I would rather no recode it. Although two different codes might not be bad anyway just to verify once it does work. The other thing is memory usage; still large. I have an idea of using a paging technique to only load bits of the mesh at a time but I really have no experience with that.

To give an idea of what it does I made a move of the simulation. The species is electrons with the source at (3, 0, 0) and there is an electric field radially outward. So the electrons should go toward the center. Instead they still seem to bounce around weirdly.

http://www.andromedaspace.com/files/simbad.mpeg

I think maybe I will implement a monte carlo integration see what that does.
Carter

Solo
Posts: 261
Joined: Fri Jun 29, 2007 12:12 pm
Location: Wisconsin

Post by Solo »

Carter: I'm interested in seeing whether I might be able to pick this up and run with it. I don't understand your grid method, and I was thinking of trying it with a plain ol' cartesian grid. I take it that this would only make the memory issues worse since the grid can't be adjusted in size. Perhaps implementing Indrek's 1/48th division could.

Also, do you have any resources I might look at?

[EDIT: how do you download the files from sourceforge? I couldn't find them there.]


It seems to me that our modeling options are constrained by the fact that we want to find the actual 3-D electric potential, (rather than assuming it) which means keeping track of the electron and ion distributions independently. I think that means we really do need a 6-D solver or 3-D PIC, both of which have large memory requirements. I think it has been suggested that the gyromotion approximation could be used to reduce the velocity to 2 components (parallel & perpendicular to the B-field) but that would get yucky since the B-field is supposed to be changed by the plasma. Then you and Art Carlson were talking about using a distribution function at the points, rather than a grid in velocity space?

The solver seems more elegant to me, but I'd bet PIC is less work.

kcdodd
Posts: 722
Joined: Tue Jun 03, 2008 3:36 am
Location: Austin, TX

Post by kcdodd »

Yes I put it under the gpl so others can work with it.

The good of a "square" grid is that you don't need to remember what the shape is so that space is saved. You could probably use an adaptive regular grid, which seems to be more common. Having sub-grids which are "square". That way you benifit from focusing on specific areas but also not so much overhead on the adaptive part. However, a lot would have to be recoded because it was done using a 6-D simplex, and not a 6-D cube. Though you could always just divide the cubes into simplices.

If you don't see the "code" tab, then hit "more". Then under "code" you can see how to access it via cvs.
Carter

Solo
Posts: 261
Joined: Fri Jun 29, 2007 12:12 pm
Location: Wisconsin

Post by Solo »

ahh, thanks, I didn't notice that set of links up there.

hmm, I guess it doesn't have to be simplices(?) to be adaptive, hadn't thought of that. I reckon it wouldn't be that much trouble to recode, seems like you did most of the work already, what with writing the solver and setting up the data handling+output (I hate I/O!).

I'll take a look at it and see what I can figure out.

[Edit: that is a big project. I think I'm out of my depth, unfortunately.]

rcain
Posts: 992
Joined: Mon Apr 14, 2008 2:43 pm
Contact:

Post by rcain »

Hi kcdodd

i like your simplex-manifold approach - personally, i would be quite interested in persevering with that approach. i have a feeling its going to be mush more flexible and powerful in the long run.

have you thought of writing the data interfaces and manipulation/transformations using geometric algebra?

(you may need to aquire a Greek keyboard ;))

btw - did you get you divide-by-zero exception to work - looked to me like it was doing what you asked it.

kcdodd
Posts: 722
Joined: Tue Jun 03, 2008 3:36 am
Location: Austin, TX

Post by kcdodd »

Yes I fixed that one, at least. No thoughts on the other stuff just wanted to get it to work first.
Carter

alexjrgreen
Posts: 815
Joined: Thu Nov 13, 2008 4:03 pm
Location: UK

Post by alexjrgreen »

Apparently chrismb is a past Research Fellow in Computational Electro-magnetics who did research on using modelling to predict sensitivity in systems and to understand when the 'theory' and a consequent electromagnetic 'model' create an unrealistic simulation of a real system.

Seems like a useful resource...
Ars artis est celare artem.

Post Reply