hayter_msa S(Q) - error returns (Trac #858)

User Paul Fitzgerald (Univ Sydney) reports that the function returns a value -1 as an error code. Should this be say NAN to avoid confusing the fitting! Example cases below. Will put in another ticket to address the instability issue.

On Tue, Feb 21, 2017 at 8:31 PM, Paul !FitzGerald ( paul.fitzgerald at sydney.edu.au ) wrote:

Hello SASView help,

I've occasionally run into an instability in the Hayter_MSA routine (both in SASView and the original Igor routine), which I usually get around by fudging the ionic strength and charge. However, this is not possible with my current set of data where we intentionally and precisely controlled the ionic strength while looking at the particle charge.

Here is an example of the problem (more examples in the attached spreadsheet):

In hayter_msa if I set: radius = 20Å Temp = 298K diel = 78 volfrac = 0.001 [salt] = 0.001 M

then vary the charge from 0 to infinity I get:

0-2 => all OK 3-32 => returns -1 for all values (error code -3) 33-115 => all OK 116-195 => returns -1 for all values (error code -3) 196-inf => all OK

The error comes from the sqcoeff(ir) function that handles rescaling. Note that this is for the Igor code - but not the compiled C version! I haven't been unable to find the HayterMSA code in the most recent version of the SASview source code so I don't know if it is the same function.

Can you tell me if (a) the theory is just not applicable to my systems and I need to try something else, or (b) if there is an implementation problem?

I did find the following in Hayter & Hensen's 1982 paper for rescaling (which I'm still trying to work through) but didn't know if it was applying to me or not: In the weak screening limit, the MSA solution of Hayter & Penfold [3] requires increasing computational precision because of the wide dynamic range of the coefficients. Under these conditions it is more efficient to use the simpler MSA solution for unscreened charged hard spheres [10] as a reference system, and to treat the weak screening as a perturbation [16]. The presence of the effective hard core allows an application of the 'optimized random phase approximation' (ORPA) [17] to express the structure factor S(K) of the screened system in terms of the structure factor So(K) of the unscreened reference system. In practice, S(K) and So(K) turn out to be almost indistinguishable for Ka < 0.5.

'''That being said, even if it turns out that the theory is not applicable, I think that there is a problem in the way that SASview deals with this problem. When the hayter_msa returns values of -1 this is an error, but SASView returns (down in the bottom left corner) "Computation completed!" instead of "error".
''' From a fitting point of view, SASView just sees a massive increase in chi square and moves in the opposite direction. This can be a problem especially when the solution is on the other side of the instability, which can be narrow. For example, if you try the above settings with vol frac = 0.01 then the instability is only from 4 to 6.

I also think that once the underlying problem should be documented in the SASView help file. We've had students run into this problem and being completely lost as to what to do.

Sorry if this creates problems.

Regards,

Paul

Migrated from http://trac.sasview.org/ticket/858

{
    "status": "closed",
    "changetime": "2017-02-27T16:11:24",
    "_ts": "2017-02-27 16:11:24.477543+00:00",
    "description": "User Paul Fitzgerald (Univ Sydney) reports that the function returns a value -1 as an error code. Should this be say NAN to avoid confusing the fitting!  Example cases below. Will put in another ticket to address the instability issue.\n \nOn Tue, Feb 21, 2017 at 8:31 PM, Paul !FitzGerald  ( paul.fitzgerald at sydney.edu.au ) wrote:\n\nHello SASView help,\n\nI've occasionally run into an instability in the Hayter_MSA routine (both in SASView and the original Igor routine), which I usually get around by fudging the ionic strength and charge.  However, this is not possible with my current set of data where we intentionally and precisely controlled the ionic strength while looking at the particle charge.   \n\nHere is an example of the problem (more examples in the attached spreadsheet):\n\nIn hayter_msa if I set: \nradius = 20\u00c5\nTemp = 298K\ndiel = 78\nvolfrac = 0.001 \n[salt] = 0.001 M\n\nthen vary the charge from 0 to infinity I get:\n\n0-2        => all OK\n3-32      => returns -1 for all values (error code -3)\n33-115  => all OK\n116-195 => returns -1 for all values (error code -3)\n196-inf  => all OK\n\nThe error comes from the sqcoeff(ir) function that handles rescaling.  Note that this is for the Igor code - but not the compiled C version!  I haven't been unable to find the HayterMSA code in the most recent version of the SASview source code so I don't know if it is the same function.\n\nCan you tell me if \n(a) the theory is just not applicable to my systems and I need to try something else, or \n(b) if there is an implementation problem?\n\nI did find the following in Hayter & Hensen's 1982 paper for rescaling (which I'm still trying to work through) but didn't know if it was applying to me or not:\nIn the weak screening limit, the MSA solution of Hayter & Penfold [3] requires increasing computational precision because of the wide dynamic range of the coefficients. Under these conditions it is more efficient to use the simpler MSA solution for unscreened charged hard spheres [10] as a reference system, and to treat the weak screening as a perturbation [16]. The presence of the effective hard core allows an application of the 'optimized random phase approximation' (ORPA) [17] to express the structure factor S(K) of the screened system in terms of the structure factor So(K) of the unscreened reference system. In practice, S(K) and So(K) turn out to be almost indistinguishable for Ka < 0.5.\n\n'''That being said, even if it turns out that the theory is not applicable, I think that there is a problem in the way that SASview deals with this problem.  When the hayter_msa returns values of -1 this is an error, but SASView returns (down in the bottom left corner) \"Computation completed!\" instead of \"error\".  \n'''\nFrom a fitting point of view, SASView just sees a massive increase in chi square and moves in the opposite direction.  This can be a problem especially when the solution is on the other side of the instability, which can be narrow.  For example, if you try the above settings with vol frac = 0.01 then the instability is only from 4 to 6.\n\nI also think that once the underlying problem should be documented in the SASView help file.  We've had students run into this problem and being completely lost as to what to do.\n\nSorry if this creates problems.  \n\nRegards,\n\nPaul\n\n",
    "reporter": "richardh",
    "cc": "",
    "resolution": "fixed",
    "workpackage": "SasView Bug Fixing",
    "time": "2017-02-22T10:51:02",
    "component": "sasmodels",
    "summary": "hayter_msa S(Q) -  error returns",
    "priority": "minor",
    "keywords": "",
    "milestone": "SasView 4.1.0",
    "owner": "pkienzle",
    "type": "defect"
}

Trac update at 2017/02/22 12:43:27:

pkienzle changed resolution from "" to "duplicate"
pkienzle changed status from "new" to "closed"

Trac update at 2017/02/22 14:38:36:

pkienzle changed resolution from "duplicate" to ""
pkienzle changed status from "closed" to "reopened"

Trac update at 2017/02/27 05:40:41: butler changed description from:

User Paul Fitzgerald (Univ Sydney) reports that the function returns a value -1 as an error code. Should this be say NAN to avoid confusing the fitting! Example cases below. Will put in another ticket to address the instability issue.

On Tue, Feb 21, 2017 at 8:31 PM, Paul FitzGerald ( paul.fitzgerald at sydney.edu.au ) wrote:

Hello SASView help,

I've occasionally run into an instability in the Hayter_MSA routine (both in SASView and the original Igor routine), which I usually get around by fudging the ionic strength and charge. However, this is not possible with my current set of data where we intentionally and precisely controlled the ionic strength while looking at the particle charge.

Here is an example of the problem (more examples in the attached spreadsheet):

In hayter_msa if I set: radius = 20Å Temp = 298K diel = 78 volfrac = 0.001 [salt] = 0.001 M

then vary the charge from 0 to infinity I get:

0-2 => all OK 3-32 => returns -1 for all values (error code -3) 33-115 => all OK 116-195 => returns -1 for all values (error code -3) 196-inf => all OK

The error comes from the sqcoeff(ir) function that handles rescaling. Note that this is for the Igor code - but not the compiled C version! I haven't been unable to find the HayterMSA code in the most recent version of the SASview source code so I don't know if it is the same function.

Can you tell me if (a) the theory is just not applicable to my systems and I need to try something else, or (b) if there is an implementation problem?

I did find the following in Hayter & Hensen's 1982 paper for rescaling (which I'm still trying to work through) but didn't know if it was applying to me or not: In the weak screening limit, the MSA solution of Hayter & Penfold [3] requires increasing computational precision because of the wide dynamic range of the coefficients. Under these conditions it is more efficient to use the simpler MSA solution for unscreened charged hard spheres [10] as a reference system, and to treat the weak screening as a perturbation [16]. The presence of the effective hard core allows an application of the 'optimized random phase approximation' (ORPA) [17] to express the structure factor S(K) of the screened system in terms of the structure factor So(K) of the unscreened reference system. In practice, S(K) and So(K) turn out to be almost indistinguishable for Ka < 0.5.

'''That being said, even if it turns out that the theory is not applicable, I think that there is a problem in the way that SASview deals with this problem. When the hayter_msa returns values of -1 this is an error, but SASView returns (down in the bottom left corner) "Computation completed!" instead of "error".
''' From a fitting point of view, SASView just sees a massive increase in chi square and moves in the opposite direction. This can be a problem especially when the solution is on the other side of the instability, which can be narrow. For example, if you try the above settings with vol frac = 0.01 then the instability is only from 4 to 6.

I also think that once the underlying problem should be documented in the SASView help file. We've had students run into this problem and being completely lost as to what to do.

Sorry if this creates problems.

Regards,

Paul

to:

User Paul Fitzgerald (Univ Sydney) reports that the function returns a value -1 as an error code. Should this be say NAN to avoid confusing the fitting! Example cases below. Will put in another ticket to address the instability issue.

On Tue, Feb 21, 2017 at 8:31 PM, Paul !FitzGerald ( paul.fitzgerald at sydney.edu.au ) wrote:

Hello SASView help,

I've occasionally run into an instability in the Hayter_MSA routine (both in SASView and the original Igor routine), which I usually get around by fudging the ionic strength and charge. However, this is not possible with my current set of data where we intentionally and precisely controlled the ionic strength while looking at the particle charge.

Here is an example of the problem (more examples in the attached spreadsheet):

In hayter_msa if I set: radius = 20Å Temp = 298K diel = 78 volfrac = 0.001 [salt] = 0.001 M

then vary the charge from 0 to infinity I get:

0-2 => all OK 3-32 => returns -1 for all values (error code -3) 33-115 => all OK 116-195 => returns -1 for all values (error code -3) 196-inf => all OK

The error comes from the sqcoeff(ir) function that handles rescaling. Note that this is for the Igor code - but not the compiled C version! I haven't been unable to find the HayterMSA code in the most recent version of the SASview source code so I don't know if it is the same function.

Can you tell me if (a) the theory is just not applicable to my systems and I need to try something else, or (b) if there is an implementation problem?

I did find the following in Hayter & Hensen's 1982 paper for rescaling (which I'm still trying to work through) but didn't know if it was applying to me or not: In the weak screening limit, the MSA solution of Hayter & Penfold [3] requires increasing computational precision because of the wide dynamic range of the coefficients. Under these conditions it is more efficient to use the simpler MSA solution for unscreened charged hard spheres [10] as a reference system, and to treat the weak screening as a perturbation [16]. The presence of the effective hard core allows an application of the 'optimized random phase approximation' (ORPA) [17] to express the structure factor S(K) of the screened system in terms of the structure factor So(K) of the unscreened reference system. In practice, S(K) and So(K) turn out to be almost indistinguishable for Ka < 0.5.

'''That being said, even if it turns out that the theory is not applicable, I think that there is a problem in the way that SASview deals with this problem. When the hayter_msa returns values of -1 this is an error, but SASView returns (down in the bottom left corner) "Computation completed!" instead of "error".
''' From a fitting point of view, SASView just sees a massive increase in chi square and moves in the opposite direction. This can be a problem especially when the solution is on the other side of the instability, which can be narrow. For example, if you try the above settings with vol frac = 0.01 then the instability is only from 4 to 6.

I also think that once the underlying problem should be documented in the SASView help file. We've had students run into this problem and being completely lost as to what to do.

Sorry if this creates problems.

Regards,

Paul

Trac update at 2017/02/27 16:11:24:

GitHub noreply@github.com commented:

In changeset a3002beb646e39f8d84ce0ad3ede19e81903c009:

#!CommitTicketReference repository="sasmodels" revision="a3002beb646e39f8d84ce0ad3ede19e81903c009"
Merge pull request SasView/sasview#260 from SasView/ticket-858

ticket-858: hayter_msa should return NAN rather than -1 for bad model parameters. Closes #142

GitHub noreply@github.com changed resolution from "" to "fixed"
GitHub noreply@github.com changed status from "reopened" to "closed"

SasView / sasmodels

hayter_msa S(Q) - error returns (Trac #858) #142