NOAA-SWPC / IPE

Ionosphere Plasmasphere Electrodynamics Model
GNU General Public License v3.0
6 stars 10 forks source link

Prompt Penetration Electric Field Validation #24

Closed twfang closed 3 years ago

twfang commented 4 years ago

We need to validate the boundary between the high latitude electric field and the dynamo solver. The default setting for crit is (15-30 degrees). A 2nd setting with (30-35 degrees) is also tested. Once IO component is included in the latest version, we need to look at the electric field directly.

twfang commented 4 years ago

The following plots are electron density (sum of all ions) at 350 km.

20130316 (Default) ne_350km_20130316_ori

20130316 (2nd setting) ne_350km_20130316_newcrit

20130317 (Default) ne_350km_20130317_ori

20130317 (2nd setting) ne_350km_20130317_newcrit

It seems like with the 2nd setting, the daytime drift might be smaller. The nighttime drifts seem to increase and last to a later time. We need to look at the drift values to confirm that.

timfullerrowell commented 4 years ago

Thanks Tzu-Wei,

Is there one or two values for each run? Could you explain a little more what 15-30 and 30-35 mean.

Thanks, Tim

On Mon, Apr 6, 2020 at 9:32 AM twfang notifications@github.com wrote:

We need to validate the boundary between the high latitude electric field and the dynamo solver. The default setting for crit is (15-30 degrees). A 2nd setting with (30-35 degrees) is also tested. Once IO component is included in the latest version, we need to look at the electric field directly.

— You are receiving this because you were assigned. Reply to this email directly, view it on GitHub https://github.com/NOAA-SWPC/IPE/issues/24, or unsubscribe https://github.com/notifications/unsubscribe-auth/AH5BFODLTQSRYDRVTBKILCDRLHYZXANCNFSM4MCL7MRQ .

twfang commented 4 years ago

These two values are used to define which potential to take. They are colatitudes. So 15-30 means 75 degree and 60 degree. The solver will use only empirical potential > 75 deg, only dynamo potential < 60 deg, and a combination of them between 60-75 mag lat.

timfullerrowell commented 4 years ago

At which latitude is the Weimer or Heelis potential taken from to apply to the dynamo solver, in each case ? Is it just one latitude or two that are applied.

On Mon, Apr 6, 2020 at 12:05 PM twfang notifications@github.com wrote:

These two values are used to define which potential to take. They are colatitudes. So 15-30 means 75 degree and 60 degree. The solver will use only empirical potential > 75 deg, only dynamo potential < 60 deg, and a combination of them between 60-75 mag lat.

— You are receiving this because you were assigned. Reply to this email directly, view it on GitHub https://github.com/NOAA-SWPC/IPE/issues/24#issuecomment-609950137, or unsubscribe https://github.com/notifications/unsubscribe-auth/AH5BFODNIGMADLNPJP6PY33RLIKVNANCNFSM4MCL7MRQ .

twfang commented 4 years ago

The Default setting use Weimer/Heelis above 75 degree while 2nd setting use Weimer/Heelis above 60 degree.

twfang commented 4 years ago

A comparison of WAM-IPE results (Weimer 2005 at high latitudes) with the empirical electric field from Manoj's model. WAM-IPE values seem to a little larger but the fluctuations look great!

Screen Shot 2020-06-08 at 11 48 24 AM

PPEF_JRO

twfang commented 4 years ago

Showing here are the lastest results from WAM-IPE v1.0 with 5 min time step. Two settings for crit1/crit 2 used are

  1. 15/30 (default)
  2. 30/35 (test setting)

20150316 (default) 20150316_wamipev1

20150316 (setting 2) 20150316_wamipe_newcrit

20150317 (default) 20150317_wamipev1

20150317 (setting 2) 20150317_wamipe_newcrit

Overall, with larger crit1/crit2 in the model, the penetration effect is significantly reduced. This also reduces the TEC during the quite time period.

twfang commented 4 years ago

I also add the electric field from setting 2 to show the comparison with WAM-IPE v1.0 and Manoj's model. You can see with different crit1/crit2, the model (black line) shows a much smaller electric field that is more comparable with the empirical model.

ppef

twfang commented 4 years ago

The document that described Barbara Emery's dynamic boundary can be found in Naomi's drive folder.

https://drive.google.com/drive/folders/0B-k7cywYbphMWUI2SmxHYXdFVlk?usp=sharing

I have looked into the code in the TIEGCM. Below are some notes from Barbara that I would like to document it. I have looked into TIEGCM a while ago and will look into this issue again.

colath.F shows crit1=theta0+5deg, and crit2=crit1+15deg advance.F calls wei05sc.F to change phid, theta0, etc from Weimer 2005. wei05sc.F calls wei05loc (inside) to set offc-4.2 deg, dskofc=0, and theta0=bndyfitr from the setboundary routine. phid and phin are also calculated in wei05loc.F which replaces my old calccloc.F routines.

It is harder to calculate theta0 from Weimer 2001, since Weimer simplified things in his 2005 code. For Heelis, theta0 is set, so you can then get crit1 and crit2 from theta0.

twfang commented 4 years ago

There are limited data from JRO. No ISR data is available at these days, only daytime JULIA data.

PPEF_Julia_data

No error bar is included in the plot, but they are about a couple of m/s for each point. I don't think the dataset is too reliable for comparing with drifts.

timfullerrowell commented 4 years ago

Yes, not very useful, and very weak dayside upward drifts.

On Mon, Aug 10, 2020 at 10:16 AM twfang notifications@github.com wrote:

There are limited data from JRO. No ISR data is available at these days, only daytime JULIA data.

[image: PPEF_Julia_data] https://user-images.githubusercontent.com/22968399/89804960-1eafcd00-daf2-11ea-8b4d-945f083d3a43.png

No error bar is included in the plot, but they are about a couple of m/s for each point. I don't think the dataset is too reliable for comparing with drifts.

— You are receiving this because you were assigned. Reply to this email directly, view it on GitHub https://github.com/NOAA-SWPC/IPE/issues/24#issuecomment-671450443, or unsubscribe https://github.com/notifications/unsubscribe-auth/AH5BFOGCC4SJLYEZJDT5LHLSAAMNTANCNFSM4MCL7MRQ .

twfang commented 4 years ago

Using the Heelis model at high latitudes and combines with different crit1/crit2 settings, here are the results compared to Manoj's empirical model similar to what has been shown above.

heelis_JRO_crit

Looks like through changing the settings for crit1/crit2, we can largely reduce the electric field at most times. It also seems to capture some of the fluctuations shown in Manoj's model (red). The changes in the electric field are a little more complex during the main phase of the storm. We need to explore to see if we can have the cross-tail potential to adopt the temporal variation of the solar wind parameters.

Some TEC comparisons for 2015/03/17 are also shown below.

Heelis + 1st setting (crit1=15, crit2=30) TEC_20150317_crit_15_30

Heelis + 2nd setting (crit1=30, crit2=35) TEC_20150317_crit_30_35

timfullerrowell commented 4 years ago

Rather than the polar cap potential, there are empirical relationships between Kp and solar wind velocity and IMF that would put back in the time dependence in Heelis. . e.g, https://doi.org/10.1002/swe.20053 Then extract the polar cap boundary from Heelis to define theta0 and crit and crit2 as in the Solomon ppt

On Thu, Aug 20, 2020 at 1:53 PM twfang notifications@github.com wrote:

Using the Heelis model at high latitudes and combines with different crit1/crit2 settings, here are the results compared to Manoj's empirical model similar to what has been shown above.

[image: heelis_JRO_crit] https://user-images.githubusercontent.com/22968399/90817610-520d0b80-e2eb-11ea-9774-60dc387e69cd.png

Looks like through changing the settings for crit1/crit2, we can largely reduce the electric field at most times. It also seems to capture some of the fluctuations shown in Manoj's model (red). The changes in the electric field are a little more complex during the main phase of the storm. We need to explore to see if we can have the cross-tail potential to adopt the temporal variation of the solar wind parameters.

Some TEC comparisons for 2015/03/17 are also shown below.

Heelis + 1st setting (crit1=15, crit2=30)

[image: TEC_20150317_crit_15_30] https://user-images.githubusercontent.com/22968399/90818834-600f5c00-e2ec-11ea-8988-09b75af7d614.png

Heelis + 2nd setting (crit1=30, crit2=35)

[image: TEC_20150317_crit_30_35] https://user-images.githubusercontent.com/22968399/90818887-71f0ff00-e2ec-11ea-9bd5-98a8ad9b55d7.png

— You are receiving this because you were assigned. Reply to this email directly, view it on GitHub https://github.com/NOAA-SWPC/IPE/issues/24#issuecomment-677870472, or unsubscribe https://github.com/notifications/unsubscribe-auth/AH5BFOFKO7QEZ4FE5A2YBELSBV5KZANCNFSM4MCL7MRQ .

twfang commented 4 years ago

Should we go with equation 1 in this paper, just take what IPE sees as solar wind parameters at each time step? Not exactly sure what the authors provided in the end.

timfullerrowell commented 4 years ago

Yes.

Do you also have the theta0 from Heelis model for the polar cap boundary? Then use the Solomon numbers, like you had before for crit1 and crit2

[image: image.png]

On Fri, Aug 21, 2020 at 1:06 PM twfang notifications@github.com wrote:

Should we go with equation 1 in this paper, just take what IPE sees as solar wind parameters at each time step? Not exactly sure what the authors provided in the end.

— You are receiving this because you were assigned. Reply to this email directly, view it on GitHub https://github.com/NOAA-SWPC/IPE/issues/24#issuecomment-678445376, or unsubscribe https://github.com/notifications/unsubscribe-auth/AH5BFOB7YBAYMXPZBCORBVLSB3ARHANCNFSM4MCL7MRQ .

twfang commented 4 years ago

Yes, then what I will do is to use solar wind parameters to calculate Kp, Kp then be used to calculate ctpoten. theta0 is then defined by (-3.80+8.48*(ctpoten*0.1875))dtr. The value will be used to modify crit1/crit2.

I will let you know how theta0 looks like throughout the 2015 event.

timfullerrowell commented 4 years ago

We might want to check the Kp value calculated by the function looks sensible by looking at a plot of the timeseries through the March storm. We might need to make sure is stays in bounds expected by the Heelis model for ctpoten as well (0 to 9).

On Fri, Aug 21, 2020 at 3:55 PM twfang notifications@github.com wrote:

Yes, then what I will do is to use solar wind parameters to calculate Kp, Kp then be used to calculate ctpoten. theta0 is then defined by (-3.80+8.48*(ctpoten*0.1875))dtr. The value will be used to modify crit1/crit2.

— You are receiving this because you were assigned. Reply to this email directly, view it on GitHub https://github.com/NOAA-SWPC/IPE/issues/24#issuecomment-678532853, or unsubscribe https://github.com/notifications/unsubscribe-auth/AH5BFOHOPZ5VDM23HGBKRNDSB3UMVANCNFSM4MCL7MRQ .

twfang commented 4 years ago

Yes, we shall check out the Kp values first. I will keep you posted.

timfullerrowell commented 4 years ago

How are you using ctpoten? Isn't Heelis driven by Kp?

On Fri, Aug 21, 2020 at 4:15 PM twfang notifications@github.com wrote:

Yes, we shall check out the Kp values first. I will keep you posted.

— You are receiving this because you were assigned. Reply to this email directly, view it on GitHub https://github.com/NOAA-SWPC/IPE/issues/24#issuecomment-678538688, or unsubscribe https://github.com/notifications/unsubscribe-auth/AH5BFOFCU7TMWE4D4RM4M4TSB3WZLANCNFSM4MCL7MRQ .

twfang commented 4 years ago

ctpoten= 15.+15.fkp+0.8fkp**2

I think this came from TIEGCM

twfang commented 4 years ago

Using the Newell equ. defined in this paper (https://agupubs.onlinelibrary.wiley.com/doi/full/10.1002/swe.20053), we are able to obtain Kp from the time-varying solar wind parameters. In the lastest setting, kp then is used to define ctpoten. ctpoten is used to determine theta0 in heelis, which is the convection reversal boundary. We do not constrain Kp to be limited to 9 since the algorithm is providing reasonable ctpoten.

image

image

image

image

image

The current method seems to provide reasonable PPEFs. After checking the SED feature and other TEC plots, the branch (heelis) will be ready to be merged into develop.

twfang commented 4 years ago

A different method that provides ctpoten using solar wind parameters was also carried out based on the paper https://agupubs.onlinelibrary.wiley.com/doi/abs/10.1029/96JA01742.

image

image

image

The results are very similar to the previous method. So, we decided to go with the Kp method at this point.

twfang commented 4 years ago

Using the same model structure (heelis model, kp using time-varying solar wind parameters), the 2003 storm is carried out. The PPEF seems to be huge, compared to the empirical model. However, it is not clear if the data is available in the empirical model to provide meaningful results. Below are some results at Jicamarca during the event.

ppef_jro_2003

The drift velocity is shown here. Large upward and downward drifts go to 200 m/s.

dft_at_JRO

TEC maps on 11/20 during the peak storm. TEC goes up to 400 TECu at where the SED feature presents.

TEC_1120

Not sure if we have anything to valid with, perhaps MIT TEC? The good thing is the model went thought the run without any problem.

timfullerrowell commented 4 years ago

This looks good.

Yes, we have MIT TEC values to compare and there are several other studies.

The vertical plasma drift does seem a bit high and too long duration. That might be the lack of shielding with the method we are using. However the very large TEC in the SED feature is very localized and could have easily been missed in observations, and it is not unlike the values observed in the Mannucci study in the Oct storm of the same year.

On Tue, Sep 15, 2020 at 1:06 PM twfang notifications@github.com wrote:

Using the same model structure (heelis model, kp using time-varying solar wind parameters), the 2003 storm is carried out. The PPEF seems to be huge, compared to the empirical model. However, it is not clear if the data is available in the empirical model to provide meaningful results. Below are some results at Jicamarca during the event.

[image: ppef_jro_2003] https://user-images.githubusercontent.com/22968399/93252801-331b6100-f753-11ea-903d-62aac38fd38d.png

The drift velocity is shown here. Large upward and downward drifts go to 200 m/s.

[image: dft_at_JRO] https://user-images.githubusercontent.com/22968399/93252838-3e6e8c80-f753-11ea-9897-db38c8c8f521.png

TEC maps on 11/20 during the peak storm. TEC goes up to 400 TECu at where the SED feature presents.

[image: TEC_1120] https://user-images.githubusercontent.com/22968399/93253416-0f0c4f80-f754-11ea-8e7e-0cd349bac1d5.png

Not sure if we have anything to valid with, perhaps MIT TEC? The good thing is the model went thought the run without any problem.

— You are receiving this because you were assigned. Reply to this email directly, view it on GitHub https://github.com/NOAA-SWPC/IPE/issues/24#issuecomment-692917938, or unsubscribe https://github.com/notifications/unsubscribe-auth/AH5BFOCETJRF5ZRIXKHN5ALSF63JPANCNFSM4MCL7MRQ .

twfang commented 4 years ago

@timfullerrowell

If you have the SED results from the 2015 and 2003 storms, could you please post them here? Also, do you have a plot from MIT TEC for 2003? Let's all put them here together. Thank you!

twfang commented 4 years ago

Some validation for 2015 and 2003 using the new Kp and the Heelis provided by Tim. https://drive.google.com/file/d/1Ga-HjtxVLNrW2EJFIU1zMQgGXZ6_a4Zo/view

twfang commented 4 years ago

Results from the lastest and the most promising method: Using Heelis at high-latitudes, but drive heelis with cross-polar cap potential (CPCP) from Weimer result at each timestep. The crit1 is theta0+5 degree and crit2 is theta0+20.

TEC on 11/20/2003

Screen Shot 2020-10-16 at 10 55 35 AM

CPCP from Weimer used in the run

Screen Shot 2020-10-16 at 10 57 48 AM

The electric field at JRO compared to PPEF empirical model on 11/19 and 11/20.

Screen Shot 2020-10-16 at 11 10 22 AM

I will test out the 2015 storm.

timfullerrowell commented 4 years ago

Looks much better. Do you have the netcdf?

On Fri, Oct 16, 2020 at 11:15 AM twfang notifications@github.com wrote:

Results from the lastest and the most promising method: Using Heelis at high-latitudes, but drive heelis with cross-polar cap potential (CPCP) from Weimer result at each timestep. The crit1 is theta0+5 degree and crit2 is theta0+20.

TEC on 11/20/2003 [image: Screen Shot 2020-10-16 at 10 55 35 AM] https://user-images.githubusercontent.com/22968399/96288433-a7d6fa80-0fa0-11eb-9207-3f90982df428.png

CPCP from Weimer used in the run [image: Screen Shot 2020-10-16 at 10 57 48 AM] https://user-images.githubusercontent.com/22968399/96288486-bf15e800-0fa0-11eb-9b6a-83bdbe7f3e5e.png

The electric field at JRO compared to PPEF empirical model on 11/19 and 11/20. [image: Screen Shot 2020-10-16 at 11 10 22 AM] https://user-images.githubusercontent.com/22968399/96288566-d7860280-0fa0-11eb-9117-56492122c665.png

I will test out the 2015 storm.

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/NOAA-SWPC/IPE/issues/24#issuecomment-710249823, or unsubscribe https://github.com/notifications/unsubscribe-auth/AH5BFOHZTJRTFHVUYVRR6L3SLB5TTANCNFSM4MCL7MRQ .

twfang commented 4 years ago

The netcdf files from this run can be found in the link below.

https://drive.google.com/drive/folders/1Hb9YIkSAe9C-oqkRHg3XZRpNTNjVzWOM?usp=sharing

timfullerrowell commented 4 years ago

Thanks I'll take a look.

On Fri, Oct 16, 2020 at 2:43 PM twfang notifications@github.com wrote:

The netcdf files from this run can be found in the link below.

https://drive.google.com/drive/folders/1Hb9YIkSAe9C-oqkRHg3XZRpNTNjVzWOM?usp=sharing

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/NOAA-SWPC/IPE/issues/24#issuecomment-710595918, or unsubscribe https://github.com/notifications/unsubscribe-auth/AH5BFOEVGHMTIDHGFLWX2JTSLCV6BANCNFSM4MCL7MRQ .

twfang commented 4 years ago

Same setting (using Weimer CPCP to drive Heelis) for 2015 storm case.

TEC on 2015/03/17

Screen Shot 2020-10-18 at 9 35 18 AM

Electric field compared to the empirical model

Screen Shot 2020-10-18 at 9 32 18 AM

I think this is good. The penetration effect seems to be reasonable. Compare to the results back in Aug, the values are not too much overestimated.

Files can be accessed in my google drive folder https://drive.google.com/drive/folders/1e--XPTK2Zb_6tzTR0BuhlsfZdPnlF2UD?usp=sharing

twfang commented 4 years ago

We now are waiting for the DMSP validation from Mariangle and Tim. If the results look good, we will merge the branch into the operational version.