Thursday, 23 August 2012

One problem found in the simplest Advection term: Differential Rotation

The last term in the Master Equation,

is the advection term in the toroidal (phi) direction. This term is dealt with the Lax-Wendroff Scheme in the complex (and more accurate) code.

Simple first order accurate Euler scheme goes haywire after a few iterations, and does not conserve the peak of magnetic field (due to numerical diffusion).

But, the Lax-Wendroff scheme, which is second order accurate improves over the Euler Scheme, and preserves the peak of magnetic field, with negligible numerical diffusion.

A problem however, generally found with Lax-Wendroff schemes is, that some ripples can be seen behind the travelling wave. The total flux, and the magnetic field are still accurately conserved. But, due to the ripples, the unsigned flux suffers change with time; which is unphysical.

 The simulation of a meridian-like flux tube shows the formation of these ripples.

For a movie of the simulation, play the video named Flux tube1.mp4 after clicking here .

Wednesday, 15 August 2012

Problem with the Euler Code

While the simulation in last post had sunspots at 22.5 degrees latitude, and their flux was not carried over 80 degrees, the time step of 100 seconds gave us quite accurate flux conservation.

But, for a sunspot placed very near the pole, we need a time step of 0.1 second for a stable numerical evolution. When run with a time step of 0.1 second, the code will approximately run for about 3 years on my laptop to complete the simulation of one solar cycle!

We seriously need to parallelise this code, if it passes all the tests we subject it to.

We also need to verify the solution it gives us, and the timespan, or the simulation time for which the results are reliable and numerically stable.

The Euler Code: Solved Imbalance..but Uncertain Accuracy

The Euler method of solving differential equations is less accurate in predicting the solution than those used in the previous code. This required us to use a smaller time step of 100 seconds for time advancing.

First, I tried only the meridional flow without any diffusion or differential rotation disturbing the magnetic field. And to my surprise, I found, that not only there was no imbalance between the hemispheres, but even the flux within each hemisphere was conserved individually!

It is logically obvious, that for a symmetric initial condition, the flux distribution at any following time should be symmetric across the equator. And so it was!

Then I gave a try to the code with only differential rotation. Still no imbalance, and accurate flux conservation!

The next step was diffusion. The diffusion only code also conserved the flux individually in hemispheres, with almost no imbalance. Sometimes, there would be an imbalance, about 10 million times less than the flux in a sunspot. But then it would become zero again.

And then, I gave a run to the code, with all the 3 effects combined, and following is a screenshot showing:
1st column: Southern Hemisphere Flux
2nd column: Northern Hemisphere Flux
3rd column: Imbalance.

I also took the liberty of making a movie of the simulation, to see the flux reach the poles in 3D! The 3D plotting was one of the best things I did while in Montana.. :)

This movie was made from a simulation that runs for about 2 years in the code, and ran for about 24 hours on my laptop!


The video quality degraded while uploading it here. For 1080p Full HD version of the video, click here .

The Imbalance Problem

Finally, some results are here, after my workstation laptop suffered 3 reinstallations of all operating systems in 3 days!

Its good to be back after an awesome summer in Montana.

A lot of new things were tried after the last post. The major concern, however was "THE IMBALANCE" of flux which we observed.

We observed that the flux in a sunspot increases/decreases rapidly as it moves around and diffuses. The problem still hasn't been completely resolved yet.

But, to see what is wrong with the code, I tried writing a simpler version of the code following the simple Euler method of solving differential equations and simple finite difference scheme.

It was thought, that this code will be quite accurate atleast for the first few iterations, and can be used to calliberate and validate the results from the full fledged code with various numerical schemes. For simplicity, I will call the simpler code (new) the 'Euler Code' and the previous code the 'Complex Code'.

Here is an example of the problem of Imbalance:
With one sunspot each in both the hemispheres, the following screenshot shows the flux values evolving over time.
1st column: Southern Hemishpere Flux
2nd column: Northern Hemisphere Flux
3rd column: Flux imbalance between the 2 hemispheres.