image357 / BoundedFromBelow

Mathematica package to check boundedness of general Higgs potentials.
BSD 3-Clause "New" or "Revised" License
2 stars 0 forks source link

Custom potential #2

Closed image357 closed 3 years ago

image357 commented 3 years ago

I added a simple example to BFB and tried to run it, At first I used your definition of polynomials

polyn and I got this error err1

and in terminal;

(base) asus@layla:~/BoundedFromBelow$ mathematica 
/usr/share/Macaulay2/Core/enginering.m2:451:15:(1):[9]: error: can't promote number to ring
/usr/share/Macaulay2/Core/enginering.m2:459:38:(1):[8]: --back trace--
/usr/share/Macaulay2/Core/classes.m2:58:53:(1):[7]: --back trace--
/usr/share/Macaulay2/Core/enginering.m2:641:40:(1):[6]: --back trace--
script.m2:3:47:(3):[5]: --back trace--
/usr/share/Macaulay2/Core/setup.m2:113:27:(1):[5]: --back trace--
/usr/share/Macaulay2/Core/setup.m2:106:18:(1):[4]: --back trace--
/usr/share/Macaulay2/Core/setup.m2:126:41:(1):[3]: --back trace--
/usr/share/Macaulay2/Core/setup.m2:129:6:(1):[2]: --back trace--
/home/dan/src/M2/M2.git/M2/Macaulay2/m2/startup.m2.in:547:63:(0):[1]: --back trace--
/home/dan/src/M2/M2.git/M2/Macaulay2/m2/startup.m2.in:668:6:(0): --back trace--
error in RunMacaulay2Script: script ended with bad exit code

Then I used my definition of polynomials; mydef

and I got this error myerr

and in terminal;

(base) asus@layla:~/BoundedFromBelow$ mathematica 
/usr/share/Macaulay2/Core/methods.m2:37:35:(1):[18]: error: no method found for applying ring to:
                              3           2
     argument   :  {{- 4eig*x1  + lsh*x1*s }} (of class List)
/usr/share/Macaulay2/Core/methods.m2:133:19:(1):[17]: --back trace--
/usr/share/Macaulay2/Core/matrix1.m2:182:75:(1):[14]: --back trace--
/usr/share/Macaulay2/Core/matrix1.m2:186:11:(1):[13]: --back trace--
/usr/share/Macaulay2/Core/matrix1.m2:281:35:(1):[12]: --back trace--
/usr/share/Macaulay2/Core/methods.m2:119:80:(1):[11]: --back trace--
/usr/share/Macaulay2/Resultants.m2:82:43:(2):[9]: --back trace--
/usr/share/Macaulay2/Core/methods.m2:119:80:(1):[7]: --back trace--
/usr/share/Macaulay2/Core/option.m2:16:8:(1):[6]: --back trace--
script.m2:4:18:(3):[5]: --back trace--
/usr/share/Macaulay2/Core/setup.m2:113:27:(1):[5]: --back trace--
/usr/share/Macaulay2/Core/setup.m2:106:18:(1):[4]: --back trace--
/usr/share/Macaulay2/Core/setup.m2:126:41:(1):[3]: --back trace--
/usr/share/Macaulay2/Core/setup.m2:129:6:(1):[2]: --back trace--
/home/dan/src/M2/M2.git/M2/Macaulay2/m2/startup.m2.in:547:63:(0):[1]: --back trace--
/home/dan/src/M2/M2.git/M2/Macaulay2/m2/startup.m2.in:668:6:(0): --back trace--
error in RunMacaulay2Script: script ended with bad exit code

Thanks Layla

Originally posted by @Layla-Kalhor in https://github.com/image357/BoundedFromBelow/issues/1#issuecomment-731523281

image357 commented 3 years ago

Before, you posted and removed this image: polyn

Here you are trying to use a rational coefficient like

lh = 0.129

within an integer ring: GetResultant[..., IntegerRing -> True]

This explains your error

/usr/share/Macaulay2/Core/enginering.m2:451:15:(1):[9]: error: can't promote number to ring

The different options that you can use are listed within the corresponding .m files in BFB.


Considering the second script: You are constructing the list of polynomials yourself. There are two problems here:

  1. You only divide by 4 where it should be 4! = 24
  2. You are not constructing a list of polynomials but a nested list of lists (compare Test.nb polynomials to yours)

The second point explains your error:

/usr/share/Macaulay2/Core/methods.m2:37:35:(1):[18]: error: no method found for applying ring to: 3 2 argument : {{- 4eigx1 + lshx1*s }} (of class List)


Two general remarks:

This will make debugging much more easy.

image357 commented 3 years ago

@Layla-Kalhor Another problem is that you are not constructing the polynomials in the right way. For instance Q5 is wrong because you don't treat the second order term in your potential correctly:

ms * s^2 / 2
Layla-Kalhor commented 3 years ago

OK, Thank you for your time,

This is a code for a simple potential;

SetDirectory[NotebookDirectory[]]
Get["BFB.m"]

H = 1/Sqrt[2] {{x1 + I*x2}, {h + I*x3}} ;    
H // MatrixForm;
HD = ConjugateTranspose[H] // ComplexExpand;
HD // MatrixForm;
HDH = FullSimplify[ HD.H] // ComplexExpand;

potential = ms*s^2/2 + ls*s^4/4 + lsh*HDH*s^2
variables = {h, x1, x2, x3, s};
Q = D[potential, {variables, 4}]/4!;
polynomials = 
 Table[Sum[
    Q[[ii, jj, kk, ll]]*variables[[jj]]*variables[[kk]]*
     variables[[ll]], {jj, 1, Length[variables]}, {kk, 1, 
     Length[variables]}, {ll, 1, Length[variables]}] - 
   eig*variables[[ii]]^3, {ii, 1, Length[variables]}]
parameters = {ms, ls, lsh, eig};

retval = GetResultant[polynomials, variables, parameters]

eigvals = 
 DeleteDuplicates[eig /. Solve[retval[["Resultant"]] == 0, eig]]

testequ = 
 Table[(polynomials /. eig -> eigvals[[4]])[[ii]] == 0, {ii, 1, 
   Length[variables]}]

Solve[testequ, variables]

And this is my error;

(base) asus@layla:~/BoundedFromBelow$ mathematica 
/usr/share/Macaulay2/Core/robust.m2:98:26:(1):[6]: error: no method for adjacent objects:
--                                                                                                   1     2                    . (of class List)
--            {{{{{{0, 0, 0, 0, 0}, {0, 0, 0, 0, 0}, {0, 0, 0, 0, 0}, {0, 0, 0, 0, 0}, {0, 0, 0, 0, --lsh*h x1}}, {{0, 0, 0, 0, .
--                                                                                                  12                          .
--    SPACE   [[1, 2, 1, 1]] (of class Array)
script.m2:3:59:(3):[5]: --back trace--
/usr/share/Macaulay2/Core/setup.m2:113:27:(1):[5]: --back trace--
/usr/share/Macaulay2/Core/setup.m2:106:18:(1):[4]: --back trace--
/usr/share/Macaulay2/Core/setup.m2:126:41:(1):[3]: --back trace--
/usr/share/Macaulay2/Core/setup.m2:129:6:(1):[2]: --back trace--
/home/dan/src/M2/M2.git/M2/Macaulay2/m2/startup.m2.in:547:63:(0):[1]: --back trace--
/home/dan/src/M2/M2.git/M2/Macaulay2/m2/startup.m2.in:668:6:(0): --back trace--
error in RunMacaulay2Script: script ended with bad exit code

Thanks for the tips, I follow them...

image357 commented 3 years ago

Again, your polynomials variable is not a list of polynomials but rather a list of lists:

/usr/share/Macaulay2/Core/robust.m2:98:26:(1):[6]: error: no method for adjacent objects: -- 1 2 . (of class List) -- {{{{{{0, 0, 0, 0, 0}, {0, 0, 0, 0, 0}, {0, 0, 0, 0, 0}, {0, 0, 0, 0, 0}, {0, 0, 0, 0, --lsh*h x1}}, {{0, 0, 0, 0, . -- 12 .

This time it's because your potential is not a scalar expression. HDH is a list containing a list containing an expression:

H = 1/Sqrt[2] {{x1 + I_x2}, {h + I_x3}} ; H // MatrixForm; HD = ConjugateTranspose[H] // ComplexExpand; HD // MatrixForm; HDH = FullSimplify[ HD.H] // ComplexExpand;


Formatting: Use triple back ticks for code blocks: ```\<optional language identifier> \<put your code here> ```

For instance ```mathematica a = 1 + 2 * I ``` to get

a = 1 + 2 * I

and single back ticks `\<code>` for inline code

Layla-Kalhor commented 3 years ago

I changed my potential; and instead of HDH, I used "(x1^2 + x2^2 + x3^2 + h^2)". I think it works, but...

SetDirectory[NotebookDirectory[]]
Get["BFB.m"]

potential = 
 ms*s^2/2 + ls*s^4/4 + lsh*s^2 * (x1^2 + x2^2 + x3^2 + h^2)
variables = {h, x1, x2, x3, s};
Q = D[potential, {variables, 4}]/4!;
polynomials = 
 Table[Sum[
    Q[[ii, jj, kk, ll]]*variables[[jj]]*variables[[kk]]*
     variables[[ll]], {jj, 1, Length[variables]}, {kk, 1, 
     Length[variables]}, {ll, 1, Length[variables]}] - 
   eig*variables[[ii]]^3, {ii, 1, Length[variables]}]
parameters = {ms, ls, lsh, eig};

retval = GetResultant[polynomials, variables, parameters]

eigvals = 
 DeleteDuplicates[eig /. Solve[retval[["Resultant"]] == 0, eig]]

testequ = 
 Table[(polynomials /. eig -> eigvals[[4]])[[ii]] == 0, {ii, 1, 
   Length[variables]}]

Solve[testequ, variables]

But it takes about an hour to run, and I had this massage; end

I have two questions,

  1. In your example, Test.nb, there are a semester for Q and a semester for polynomials. If I want to add new potential, these two semesters remain the same and I should not change, right?
    
    Q = D[potential, {variables, 4}]/4!;

polynomials = Table[Sum[ Q[[ii, jj, kk, ll]]variables[[jj]]variables[[kk]] variables[[ll]], {jj, 1, Length[variables]}, {kk, 1, Length[variables]}, {ll, 1, Length[variables]}] - eigvariables[[ii]]^3, {ii, 1, Length[variables]}]



2. In the last line of the Test.nb, Can the equation be solved based on parameters? And get the range of parameters to bound the potential from below?
instead of "Solve[testequ, variables]" use "Solve[testequ, parameters]"?

I want to know how the parameters change so that the potential is bounded from below?

Best Wishes
Layla
image357 commented 3 years ago

General answer

Have a look at the publications and try to understand the algorithm that governs the boundedness test. Especially, the second half of my thesis. You will then understand the purpose of the expressions in Test.nb. I will try to summarize it here:

  1. Derive the fourth rank tensor Q.
  2. Calculate eigenvalues via the resultant method.
  3. Potential is bounded iff for all real eigenvalues, all non-trivial, real solutions (i.e. eigenvectors) belong to a positive eigenvalue.

The last point gives you a constraint on your parameters.

Speaking differently: The potential is not bounded from below iff a real and non-positive (i.e. zero or negative) eigenvalue has a real solution different than x=y=...=0.

Specific answers

But it takes about an hour to run, and I had this massage;

That's normal for symbolic evaluations. The resultant algorithms implemented in Macaulay2 are not efficient (think about calculating the determinant of a huge matrix symbolically and not numerically).

In your example, Test.nb, there are a semester for Q and a semester for polynomials. If I want to add new potential, these two semesters remain the same and I should not change, right?

What do you mean by semester?

In the last line of the Test.nb, Can the equation be solved based on parameters? And get the range of parameters to bound the potential from below? instead of "Solve[testequ, variables]" use "Solve[testequ, parameters]"?

See the general answer from above. Once you have Solve[testequ, variables] you can test for positivity by checking for real solutions and positive eigenvalues (evaluate solutions and eigenvalues with a specific parameter set).

I want to know how the parameters change so that the potential is bounded from below?

You can scan parameters for positivity with PositivityTest. It supports a semi-numerical approach. All parameters must be numbers. However, it will treat them as rationals. This is faster than a symbolic calculation but still not as fast as a full numerical treatment would guarantee. This algorithmic caveat is based on the possible evaluations within Macaulay2 (to my knowledge the only computer algebra system that supports resultants for multivariate polynomials).

Summary

The runtime you experience is expected. Mathematica wants to tell you that something is taking very long. You can use PositivityTest to check parameters or follow points 1 - 3 (which PositivityTest does, too). Currently, the algorithm scales very badly with the number of variables. A full numerical treatment would weaken that caveat. The algorithms for that are described in my publications and references therein, but the effort to implement them should not be underestimated.

Lastly, for general higgs potentials (of many variables), one can prove that exact positivity cannot be decided based on "simple" inequalities (contrary to the case for most 2HDMs).

Layla-Kalhor commented 3 years ago

Thank you so much,

In your example, Test.nb, there are a semester for Q and a semester for polynomials. If I want to add new potential, these two semesters remain the same and I should not change, right?

What do you mean by semester?

I mean, these two definitions (Q and polynomials) should not be changed for other models, I just have to change the potential.

Layla-Kalhor commented 3 years ago

I thank you and Prof. Dr. M. Margarete Mühlleitner, and her perfect group. I thank Dr. Igor P. Ivanov for introducing you to me.

Your information is really useful, and thank you for answering my questions patiently. However, I will ask the question again!!

Prof. Dr. M. Margarete Mühlleitner has great students such as you and Phillip Basler, ... . We have used the BSMPT package in our works and we have new ideas about it. The authors of this article, like you, patiently answer our questions. I hope the best for your group...

image357 commented 3 years ago

I mean, these two definitions (Q and polynomials) should not be changed for other models, I just have to change the potential.

Yes, that's correct.

Your information is really useful, and thank you for answering my questions patiently. However, I will ask the question again!!

Prof. Dr. M. Margarete Mühlleitner has great students such as you and Phillip Basler, ... .

You are welcome. I will let them know ;) I'll close this issue now but feel free to open a new one if you run into any problem. Cheers, Marcel


PS: Here is some example code for PositivityTest

PositivityTest[2*x^4 - 1*x^2*y^2 + 2*y^4, {x, y}]
PositivityTest[-2*x^4 + 1*x^2*y^2 + 2*y^4, {x, y}]
PositivityTest[0*x^4 + 1*x^2*y^2 + 2*y^4, {x, y}]
Layla-Kalhor commented 3 years ago

Thanks