February 2, 2011 (2/2/2011) Begin modifications to the zaphod version of smeiipstd0n_intel To run type ./smeiipstd0n_intel smei=/smeitom/int/inter_data/base/*.txt nagoya=nagoya,,yearly The mk_smeiipstdn_intel program just compiled this program version. Note - this program works if the shell has an unlimited ability to read files. When the directory has more than ~32,768 entries (points files), this requires "ulimit -s unlimited" to be set in unix?, maybe some other command in another operating system. On 2/2/2011 John Clover discovered that the mkpostdn.f subroutine was not the latest version that includes the ability to correctly handle leap years. The most current version is now included in the subroutine set that will now be used to provide the we-based tomography. On 2/3/2011 on tests, it turns out that the current tomographic run on CR 2092 now begins a day earlier. Thus, this fix is imporntant for the times when the start of the tomography is just at the beginning of the year for non leap year rotations. Thus, the correct version of mkpostdn.f fixes a problem present when at the beginning of 2010 - namely for the January 17 2003 CME where the small density enhancement now crosses the STEREO-B location centered at 01/20/2010 20:00UT as opposed to 01/19/2010 00:00UT (with a 400kms default speed) using the old wrapper-produced set of images. On 2/4/2011 I have set different limits on a couple parameters to keep the inner portion of the volume to be "eaten" at the end of the run when archival data are being used. To do this, I first set NTmaxx to 60 (It was 55 days). I then tested this to see if this made a difference. It didn't; The "eaten away" part was identical when I did this for CR2099. I returned the value of NTmaxx to 55. I then set the value of NLmaxG to 400,000, and the beginning orbit: ibtt(1) = ((XCbeg + NCoff) - XCORBEG)*27.2753 ! Approx. days since the beginning orbit (intel >02/02/11 BVJ) and the ending orbit: iett(1) = ((XCend + NCoff) - XCORBEG)*27.2753 ! Approx. days since the ending orbit (intel >02/02/11 BVJ) The former was (XCbeg + 0.15 + NCoff) The latter was (XCbeg + NCoff + 1.3 + xInc) This new change allows a 56 more orbits to be read at the rotation beginning and a total of 780 orbits to be read. The latter was the main reason the current program using archival data was "eaten" at the end, and this new fix seems to allow a continuous blend of one rotation to another for the current (default) speeds set in the tomography program. Thus, on 2/4/2011 I have made the main changes to: smeiipstd0n_intel.f as I stated above. I also made small changes to: On 2/4/2011 I smeiread.f to print out a more concise statement of the read and accumulated points files, and I also made small changes to: read_smei8.f in order to provide a statement at the time location xTT where the number of files read exceeds the MJDCar total limit allowed. The value of MJDCar is set to MJDref = JDCar(I)-Sun__MJDtoJD+(XCtest2-I)*(JDCar(I+1)-JDCar(I)) where I = XCtest2, the value of the last Carrington Rotation to be read in the main program. This translates to a modified Julian date time of the value of XCTest2, which is the last Julian Date allowed to be read minus a small time. To get this value, the Julian Date of the last time to be read (JDCar(I) of the value of XCTEST2) is converted to Modified Julian data by subtracting the constant Sun__MJDtoJD = 2400000.50000000 (double precision) and adding a zero (XCtest2-I). Humm, at one time I expect XCtest2 did not equal I. In any case this value is now checked against a value of the time in Modified Julian Date read from the points file input in subroutine read_smei8.f. Evidently, at one time there were bogus time values present in the input that were eliminated by this check. After this check the program continues to read files up to the MX limit. In the current scheme for CR2099, no points files read exceed the xTT check, or beyond orbit 396566 (the last orbit now read is 396561). The first time the xTT limit is exceeded, if it is exceeded, the program now writes this out (in read_smei8.f) to flag that there may be some error in the program. For CR2099 there are 329877 accumulated points (lines-of-sight) read and kept in smeiread.f, but only 261876 accepted through the first set of criteria in read_smei8.f, and transferred to the main program smeiipstd0n_intel for further checks and analysis. On 2/4/2011 I also tried the current program compiled on zaphod to see if the rotation 2025 was correct. This rotation begins on January 2 following a leap year, and at one time was displaced by a day from what it should have been for the event on January 20, 2005. This event has more than 400,000 points files accepted. I have now set the program limit to 500,000. I noticed that the print statement in read_smei8.f did not work to tell users that the program had reached its limit. I also noticed that the subroutine did not stop writing into memory when this happened for the last few LOSs of the points files. An if statement in this subroutine now hopefully fixes this problem so that the program does not write into protected memory, and this fix now tests OK. The original test to check that the leap year is handled OK, now also shows that this is the case, or at least for the January 20, 2005 CME shock that arrives at Earth on January 21 at ~18 UT the analysis gives a value that is nearly identical to that on the Web that shows the CME shock density enhancement arriving at 18 UT on the 21st. If the Web analysis was done correctly by making two runs, one that began before the beginning of the year, then the new analysis is indeed OK. On 2/6/2011 I ran the smeiipstd0n 32-bit program compiled on new.cass183.ucsd.edu and tested it run on zaphod with the same points files as the smeiipstd0n_intel program compiled on zaphod. I did this for CR 2099, and checked the program using the IDL-produced ecliptic cuts. The program compiled on new.cass183, and the program compiled on zaphod show a few minuscule differences, but are otherwise so close to identical, that there is no reason not to use the zaphod-compiled tomography program. At the end of the program, the inner portion of the volume was "eaten" even more than the earlier compilations on zaphod. There is little wonder about this, since the smeiipstd0n program on cass183 can accept only fewer lines of sight than can this early compilation. (2/7/2011) smeiipstd0n_hr On 2/7/2011 I produced another version of smeiipstd0n_intel.f called smeiipstd0n_hr with its own mkfile called mk_smeiipstd0n_hr. I first upped the value of NGREST to 4 (from 2) and upped the value of NGRESS to 6 (from 3). I will now see if this higher resolution version of the program runs (without IPS data), how much space it takes on zaphod, and if it gives sensible results. The program began at about 5:30 am and it seems to take a value of about 1819172 memory allocation. On 2/8/2011 I found this program does not seem to run with NGREST different from NVREST. I found several errors when this was done, the major one was that the values of velocity errors were not properly placed into the routines used to produce the errors or the routines to provide the nv3 files at high resolution. The arrays that went into these files were not dimensioned properly. When they are dimensioned properly, the program still did not work to provide high resolution files. An error was produced by the program that wrote the error file. An error was also produced in an attempted write of the high resolution files. I am not sure why this is so, but the indication is that the copyvtovdn program is making an error when the temporal resolutions are not identical for both the density and the velocity files. On 2/9/2011 I discovered why an error was produced by the program that wrote the error file. The error file for velocity was not being intitalized before being input to the write program. Because the file was not initalized in this version of fortran, garbage remained in the file that was then transferred to the write program through the copyvtodvn program. In the past, the non-initalization produced a zero which was then flagged bad in the main program. To fix I initalized the error program with a 1.0. I suppose better would be to initalize the program with a zero, but in transfer, the zero may be changed to some low number (or the bad may be changed to something else which then may not give erronous results in the write program. In any case, this initalization has fixed the problem with both the error write and the nv3h* write, and the program now works and gives reasonable results when run. I checked the program at the NGREST = 2 and NVREST = 2 values with previous runs. The program now runs at higher resolutions with NGREST different from NVREST. On 2/9/2011 I discovered a potential error in that the higher temporal resolutions seem to take a much longer time for structures to propagate outward than when the lower resolution programs are run. Preliminary indications are that this is caused by an error in the mkshiftn0.f subroutine and attributed to the constant AKMPAUPI set up at the beginning of the subroutine. If so, this is a major error in the program, and will take significant consideraton to fix. This constant propagates through the subroutine but mostly affects the level shifts in time from one level to another, and the solar rotation amount from one level to another (the constant ERot__Fract which it modifies). I have renamed mkshiftn0.f to mkshiftn0n.f and in the main program to distinguish this version from the previous. I also added the argument "speed" in the argument list of the mkshiftn0.f (now mkshiftn0n.f) subroutine, since this was previously used as set at the beginning of the mkshiftn0.f subroutine to provide a default traceback speed (at 425 kms) rather than the speed set on input to the main program. This default is only evoked if no velocity data exists at a level location in the mkshiftn0.f (now mkshiftn0n.f) subroutine. On 2/10/2011 I found that in replacing the value of AKMPAUPI with AKMPAUPD, I was able to get the current structures in the solar wind at NGREST = 4 to move outward over a period of 4 days from the Sun for the CR2099 (August 5, 2010) analysis with an input speed of 400 kms. On 2/10/2011 I am now re-running both CR2099 and CR2092 with NGREST = 2 to see how the tomography is affected with a "normal" run of the current version of smeiipstd0n_hr. There will be considerable differences. I have also provided another input NinterDD that allows an automatic input of the preferred value of NinterD (intermediate interpolated temporal values). In the original program with NGREST = 2 and NinterD = 3 high resolution nv3h* and other files were produced every 3 hours. This is still the case, but now the program asks what you if you want other than the preferred default (not necessairly NinterD = 3) if NGRESS > 2 in order to provide 3-D temporal reconstructions at a given hourly cadence. On 2/11/2011 I looked at several cases where mkshiftn0n.f is now used. Yes there are considerable differences. As expected, the new main program now provides a solar wind model with NGREST = 2 that is twice as fast as before. For the rotation 2092, this means that the small little bump no longer arrives at STEREO-B on 21 January, but late 18 January. The identification of this as the January 14 CME is thus incorrect. The main CME observed IS the 14 January CME. The 17 January CME comes later and is what is measured as the later event at STEREO-B on 23 January. Another test uses the 1 August 2010 CME. IN this case the event arrives at Earth at the exact same time as before, and simply looks more spread out in radial extent. The other test was for Carrington Rotation 2025 (the January 20 CME) Here the event arrives at Earth at the exact same time as before, but now the run using intel2 has the shock density enhancement looking a lot like that published in Solar Wind 12 where the speed was changed to 1000kms. This is exactly as I would expect. On 2/13/2011 I checked to see that the version of the program compiled on zaphod gives the same result when run on bender. For CR2025 the results viewed by ecliptic cuts are identical. For CR2099 the program bombs in write3d_info3dM_HR.f when smeiipstd0n_hr is run using NGREST = 4 and NGRESS = 6. On 2/14/2011 I checked the three versions of the program compiled on bender (using a fortran 02, 03, and 04 option) with CR 2025 and between the compilation of smeiipstdo0n_intel2 and the first EA_2025.00 file: 02 case 52 min 36 sec 03 case 43 min 10 sec 04 case 53 min 41 sec I doubt the timing on these some, because in the first case I was not too strict about the beginning of the run relative to the compilation, and on the final case (04), the run was operating simultaneously with the 03 case. John says that the compiled version using my mk_file is not the optimal run since with his Makefile all subroutines are optimized including those in the libraries. On 2/14/2011 I learned to use the Makefile technique type, for instance, "make smeiipstd0n_hr", and this will run the current instructions in Makefile with the smeiipstd0n_hr option. If smeiipstd0n_hr has been run earlier then better delete this binary file to be sure the program is recompiled.