Full text
2013 59 Jorge Monforte García Study of the critical and low temperature properties of finite dimensional spin glasses Departamento Director/es Física Teórica Tarancón Lafita, Alfonso Ruiz Lorenzo, Juan Jesús Director/es Tesis Doctoral Autor Repositorio de la Universidad de Zaragoza – Zaguan http://zaguan.unizar.es UNIVERSIDAD DE ZARAGOZA
Departamento Director/es Jorge Monforte García STUDY OF THE CRITICAL AND LOW TEMPERATURE PROPERTIES OF FINITE DIMENSIONAL SPIN GLASSES Director/es Física Teórica Tarancón Lafita, Alfonso Ruiz Lorenzo, Juan Jesús Tesis Doctoral Autor 2013 Repositorio de la Universidad de Zaragoza – Zaguan http://zaguan.unizar.es UNIVERSIDAD DE ZARAGOZA
Departamento Director/es Director/es Tesis Doctoral Autor Repositorio de la Universidad de Zaragoza – Zaguan http://zaguan.unizar.es UNIVERSIDAD DE ZARAGOZA
Study of the critical and low temperature properties of finite dimensional spin glasses Memoria de tesis doctoral presentada por Jorge Monforte Garc ´ ıa Directores Alfonso Taranc´ on Lafita Juan Jes´ us Ruiz Lorenzo Departamento de F´ısica Te´orica, Facultad de Ciencias Universidad de Zaragoza Instituto de Biocomputaci´on y F´ısica de Sistemas Complejos Zaragoza, 29 de abril de 2013
D. Alfonso Taranc´on Lafita, Catedr´atico en el Departamento de F´ısica Te´orica de la Universidad de Zaragoza y D. Juan Jes´us Ruiz Lorenzo, Profesor Titular del ´area de F´ısica Te´orica del Departamento de F´ısica de la Universidad de Extremadura HACEN CONSTAR que la presente memoria de Tesis Doctoral presentada por Jorge Monforte Garc´ıa y titulada “Study of the critical and low temperature properties of finite dimensional spin glasses” ha sido realizada bajo su direcci´on en el Departamento de F´ısica Te´orica de la Universidad de Zaragoza y en el Instituto de Biocomputaci´on y F´ısica de Sistemas Complejos, BIFI. El trabajo recogido en esta memoria se corresponde con lo planteado en el proyecto de tesis doctoral aprobado en su d´ıa por el ´organo responsable del programa de doctorado. Y para que conste, cumplimiento de la legislaci´on vigente, informan favorablemente sobre la referida Tesis Doctoral y autorizan su presentaci´on para su admisi´on a tr´amite. Zaragoza, a 29 de abril de 2013 Los directores de la Tesis Fdo: Alfonso Taranc´on Lafita Fdo:Juan Jes´us Ruiz Lorenzo
i
A mi familia iii
iv
CONTENTS xi investigador. Agradecer a su vez a Jes´us Clemente Gallardo, Pierpaolo Bruscolini y Jos´e Luis Alonso por ense˜narme y guiarme en las que han sido las primeras clases que he impartido en la universidad. Gracias tambi´en a Joaqu´ın Sanz, por compartir tan buenos momentos prepar´andolas. En general y para concluir, quiero agradecer a todas aquellas personas que han vivido conmigo la realizaci´on de esta tesis doctoral. Desde lo m´as profundo de mi coraz´on os agradezco todo el apoyo, ´animo y colaboraci´on que me hab´eis dado, pero sobre todo, todo vuestro cari˜no, amor y amistad. Jorge Monforte Garc´ıa ha realizado esta Tesis Doctoral gracias a una beca del BIFI, una beca Predoctoral de Formaci´on de Personal Investigador del Gobierno de Arag´on, a los proyectos “Complejidad en materiales y fen´omenos de transporte” FIS2006-08533-C03-02 y FIS2009-12648-C03-02 del MICINN, “Grupo de excelencia de biocomputaci´on” del Gobierno de Arag´on E24/3 y “SCC Computing” del VII proyecto Marco de la UE.
xii CONTENTS
Resumen Esta memoria ha sido dedicada al estudio de modelos de vidrios de esp´ın con interacciones a corto alcance, en concreto el modelo de Potts de vidrios de esp´ın y el de Edwards-Anderson. El objetivo principal de esta Tesis Doctoral ha sido el estudio de las transiciones de fase que estos modelos presentan as´ı como la caracterizaci´on de su fase de vidrio de esp´ın a bajas temperaturas. La complejidad que presentan los vidrios de esp´ın exigen el desarrollo de sofisticadas herramientas para su estudio, las cuales pueden ser aplicadas en otras ramas de la ciencia como el plegamiento de prote´ınas. Para el desarrollo de esta tesis se han utilizado programas propios escritos en lenguaje C y la m´aquina dedicada Janus del BIFI, as´ı como en menor medida otras infraestructuras como el Cluster del BIFI. Un vidrio de esp´ın es una colecci´on de momentos magn´eticos, espines, que a baja temperatura presenta un estado congelado desordenado, la fase de vidrio de esp´ın. En esta fase, el sistema posee caracter´ısticas muy interesantes. Los tiempos de relajaci´on son extremadamente largos debido a un paisaje de energ´ıa muy complicado. Una de las principales causas de ello es la frustraci´on, que consiste en que los espines no son capaces de encontrar un estado estable debido a que hay competencia entre distintas interacciones con los espines vecinos. Los primeros vidrios de esp´ın que se estudiaron, en la d´ecada de 1970, fueron los vidrios de esp´ın met´alicos o can´onicos, compuestos por una base met´alica en la que se a˜naden impurezas magn´eticas. Desde entonces se han dedicado muchos trabajos al estudio tanto experimental como te´orico de los vidrios de esp´ın, aunque a´un quedan muchas inc´ognitas abiertas. Esta tesis intenta realizar una peque˜na aportaci´on a este vasto campo de investigaci´on. El Cap´ıtulo 1 es una introducci´on a los vidrios de esp´ın donde se explican qu´e son estos materiales, sus caracter´ısticas y se dan algunos ejemplos de materiales reales que presentan el comportamiento de un vidrio de esp´ın. Adem´as, se presentan varios modelos realistas de vidrios de esp´ın as´ı como aproximaciones que tienen soluci´on anal´ıticas que nos permiten lanzar hip´otesis sobre el comportamiento de los modelos m´as realistas. xiii
xiv CONTENTS En el Cap´ıtulo 2 se estudia el comportamiento cr´ıtico del modelo de Potts de vidrios de esp´ın en tres dimensiones con 5 y 6 estados. Este modelo presenta un diagrama de fases muy rico y por eso recibe bastante atenci´on. En este caso, nosotros caracterizamos la transici´on a la fase de vidrio de esp´ın y estudiamos su dependencia con el n´umero de estados, as´ı como buscamos la posible existencia de otra transici´on de fase, esta vez a una ferromagn´etica. Uno de los objetivos principales del Cap´ıtulo 3 es el estudio de ciertas caracter´ısticas de la fase de vidrio de esp´ın de sistemas finitos con interacciones de rango finito, como por ejemplo estabilidad estoc´astica, Replica Equivalence (equivalencia de r´eplicas), Overlap Equivalence (equivalencia de overlap) y ultrametricidad. Para ello se estudian las fluctuaciones entre muestras del modelo de Edwards-Anderson en tres dimensiones. En el Cap´ıtulo 4 se investiga, utilizando t´ecnicas fuera del equilibrio, si existe transici´on de fase en un vidrio de esp´ın en tres dimensiones en presencia de un campo magn´etico externo, ya que los dos principales escenarios te´oricos predicen comportamientos antag´onicos. En el Cap´ıtulo 5 aplicamos una t´ecnica alternativa para estudiar transiciones de fase en vidrios de esp´ın, el an´alisis de las singularidades complejas de la funci´on de partici´on. Esta t´ecnica, fue desarrollada en 1952 por Lee y Yand y desde entonces se ha aplicado a multitud de sistemas f´ısicos, por lo que queremos estudiar su aplicaci´on a vidrios de esp´ın. En el Cap´ıtulo 6 se estudian los fen´omenos de rejuvenecimiento y memoria que presentan los vidrios de esp´ın cuando son sometidos a cambios de temperaturas en su fase de vidrio de esp´ın fuera del equilibrio. Se intentar´a reproducir el impresionante experimento dip en el que estos fen´omenos se evidencian claramente. En el Cap´ıtulo 7 se recoge un resumen de los trabajos de investigaci´on en los que he trabajados dentro de la Janus Collaboration pero que no forman la parte principal de investigaci´on de esta Tesis Doctoral. Finalmente, el Cap´ıtulo 8 est´a dedicado a las conclusiones. Como resultado de todo este trabajo han sido realizadas las siguientes publicaciones •R. Alvarez Ba˜nos, A. Cruz, L. A. Fernandez, A. Gordillo-Guerrero, J. M. Gil-Narvion, M. Guidetti, A. Maiorano, F. Mantovani, E. Marinari, V. Martin-Mayor, J. Monforte-Garcia, A. Mu˜noz Sudupe, D. Navarro, G. Parisi, S. Perez-Gaviro, J. J. Ruiz-Lorenzo, B. Seoane, S. F. Schifano, A. Tarancon, R. Tripiccione and D. Yllanes, J. Stat. Mech. P05002 (2010). “Critical Behavior of Three-Dimensional Disordered Potts Models with Many States”. •R. A. Ba˜nos, A. Cruz, L. A. Fernandez, J. M. Gil-Narvion, A. Gordillo-
CONTENTS xv Guerrero, M. Guidetti, D. I˜niguez, A. Maiorano, F. Mantovani, E. Marinari, V. Martin-Mayor, J. Monforte-Garcia, A. Mu˜noz-Sudupe, D. Navarro, G. Parisi, S. Perez-Gaviro, F. Ricci-Tersenghi, J. J. RuizLorenzo, S. F. Schifano, B. Seoane, A. Taranc´on, R. Tripiccione and D. Yllanes, Phys. Rev. B 84, 174209 (2011). “Sample-to-sample fluctuations of the overlap distributions in the three-dimensional EdwardsAnderson spin glass”. •R. A. Ba˜nos, J. M. Gil-Narvion, J. Monforte-Garcia, J. J. Ruiz-Lorenzo and D. Yllanes, J. Stat. Mech. P02031 (2013). “Numerical Study of the Overlap Lee-Yang Singularities in the Three-Dimensional EdwardsAnderson Model”. •F. Belletti, A. Cruz, L. A. Fernandez, A. Gordillo-Guerrero, M. Guidetti, A. Maiorano, F. Mantovani, E. Marinari, V. Martin-Mayor, J. Monforte, A. Mu˜noz Sudupe, D. Navarro, G. Parisi, S. Perez-Gaviro, J. J. Ruiz-Lorenzo, S. F. Schifano, D. Sciretti, A. Tarancon, R. Tripiccione and D. Yllanes, J. Stat. Phys. 135, 1121 (2009). Eprint: arXiv:0811.2864. “An in-depth look at the microscopic dynamics of Ising spin glasses at fixed temperature”. •R. ´ Alvarez Ba˜nos, A. Cruz, L. A. Fernandez, J. M. Gil-Narvion, A. Gordillo-Guerrero, M. Guidetti, A. Maiorano, F. Mantovani, E. Marinari, V. Martin-Mayor, J. Monforte-Garcia, A. Mu˜noz Sudupe, D. Navarro, G. Parisi, S. Perez-Gaviro, J. J. Ruiz-Lorenzo, S. F. Schifano, B. Seoane, A. Tarancon, R. Tripiccione and D. Yllanes, J. Stat. Mech. P06026 (2010). Eprint: arXiv:1003.2569. “Nature of the spinglass phase at experimental length scales”. •R. ´ Alvarez Ba˜nos, A. Cruz, L. A. Fernandez, J. M. Gil-Narvion, A. Gordillo-Guerrero, M. Guidetti, A. Maiorano, F. Mantovani, E. Marinari, V. Martin-Mayor, J. Monforte-Garcia, A. Mu˜noz Sudupe, D. Navarro, G. Parisi, S. Perez-Gaviro, J. J. Ruiz-Lorenzo, S. F. Schifano, B. Seoane, A. Tarancon, R. Tripiccione and D. Yllanes, Phys. Rev. Lett. 105, 177202 (2010). Eprint: arXiv:1003.2943. “Static versus dynamic heterogeneities in the D= 3 Edwards-Anderson-Ising spin glass”. •R. ´ Alvarez Ba˜nos, A. Cruz, L. A. Fernandez, J. M. Gil-Narvion, A. Gordillo-Guerrero, M. Guidetti, D. I˜niguez, A. Maiorano, E. Marinari, V. Martin-Mayor, J. Monforte-Garcia, A. Mu˜noz Sudupe, D. Navarro, G. Parisi, S. Perez-Gaviro, J. J. Ruiz-Lorenzo, S. F. Schifano, B. Seoane, A. Tarancon, P. Tellez, R. Tripiccione and D. Yllanes,
xvi CONTENTS PNAS 109 6452 (2012). Eprint: arxiv:1202.5593. “Thermodynamic glass transition in a spin glass without time-reversal symmetry”. •M. Baity-Jesi, R. A. Ba˜nos, A. Cruz, L. A. Fernandez, J. M. GilNarvion, A. Gordillo-Guerrero, M. Guidetti, D. I˜niguez, A. Maiorano, F. Mantovani, E. Marinari, V. Martin-Mayor, J. Monforte-Garcia, A. Mu˜noz Sudupe, D. Navarro, G. Parisi, M. Pivanti, S. Perez-Gaviro, F. Ricci-Tersenghi, J. J. Ruiz-Lorenzo, S. F. Schifano, B. Seoane, A. Tarancon, P. Tellez, R. Tripiccione and D. Yllanes, Eur. Phys. J. Special Topics 210, 33 (2012). “Reconfigurable computing for Monte Carlo simulations: Results and prospects of the Janus project”.
CONTENTS xvii Nota filol´ogica El resumen, el primer cap´ıtulo y el cap´ıtulo de las conclusiones han sido escritos en espa˜nol. El resto de la tesis, incluidas las traducciones del primer cap´ıtulo y del cap´ıtulo de las conclusiones, ha sido escrito en ingl´es con el objetivo de permitir su lectura a un grupo m´as amplio de personas.
xviii CONTENTS
Cap´ıtulo 1 Introducci´on Los vidrios de esp´ın1son sistemas magn´eticos (es decir, una colecci´on de espines) que presentan una transici´on de fase a una fase congelada de baja temperatura desde una paramagn´etica. Sin embargo, esta fase no exhibe orden de largo alcance (mientras que los materiales ferromagn´eticos y antiferromagn´eticos s´ı), por lo que esta fase es una especie de desorden congelado. Por tanto, la magnetizaci´on local mies no nula mientras que la magnetizaci´on media M=∑imi N(1.1) donde Nes el n´umero total de espines, y la magnetizaci´on a momento k Mk=∑ i e−ik·rimi N(1.2) se anulan para todos los momentos k. La ausencia de un orden de largo alcance (a diferencia de los materiales antiferromagn´eticos) se puede comprobar con experimentos de scattering de neutrones. Los vidrios de esp´ın met´alicos o can´onicos fueron el primer tipo de vidrios de esp´ın estudiado. Estos materiales son aleaciones met´alicas creadas a˜nadiendo impurezas magn´eticas a una base met´alica, por ejemplo, CuMn. Es bien conocido que en un ferromagneto (como Fe), la interacci´on magn´etica es calculada con la interacci´on de canje, por lo que se obtiene que H=−JS1S2(1.3) donde S1yS2son los espines (es decir, los momentos magn´eticos) de los ´atomos magn´eticos. Sin embargo, en un vidrio de esp´ın met´alico, los ´atomos 1Se han publicado muchas reviews, nos centraremos en las Refs. [2, 3, 5]. 1
2CAP´ ITULO 1. INTRODUCCI ´ ON magn´eticos son impurezas por lo que se tiene una especie de interacci´on de canje indirecta, una impureza magn´etica interacciona con un electr´on de conducci´on quien, despu´es, interacciona con otra impureza magn´etica. Esta interacci´on es la llamada interacci´on RKKY (debido a que fue estudiada por Runderman y Kittel en 1954 [6], Kasuya en 1956 [7] y Yosida en 1957 [8]) y la expresi´on de la interacci´on es J(r) = J0 cos(2kFr+φ0) (kFr)3(1.4) donde J0yφ0son constantes y kFes el n´umero de ondas de Fermi del metal anfitri´on (en nuestro ejemplo, Cu). Los tiempos de relajaci´on en la fase congelada de vidrio de esp´ın son extremadamente largos, por lo que el estudio de la din´amica fuera del equilibrio es muy ´util para comparar con experimentos. Adem´as, algunos fen´omenos t´ıpicos de los vidrios de esp´ın aparecen en este r´egimen. El comportamiento del sistema depende del proceso de enfriamiento y el tiempo, tw, que el sistema haya estado en la fase de vidrio de esp´ın, es decir, los vidrios de esp´ın exhiben aging (ver, por ejemplo, Ref. [9]). Si el sistema evoluciona un tiempo twen una temperatura fija, T, en la fase de vidrio de esp´ın, aparecen dos ejemplos de fen´omenos de aging: la magnetizaci´on termorremanente (el sistema evoluciona en presencia de un campo magn´etico externo que luego se retira) y la magnetizaci´on de un enfriamiento a campo cero, ZFC (el campo magn´etico externo se enciende tras haber pasado el sistema un tiempo twen la fase de vidrio de esp´ın). Un par de ejemplos de este tipo de experimentos se muestra en las Figuras 1.1 y 1.2. Adem´as, si la temperatura en la fase de vidrio de esp´ın no se mantiene constante aparece otros fen´omenos, como el rejuvenecimiento y la memoria (ver Cap´ıtulo 6). En resumen, las principales caracter´ısticas de un vidrio de esp´ın en su fase de vidrio de esp´ın son que los momento magn´eticos est´an congelados, ausencia de orden de largo alcance (Mk= 0 y M= 0), tiempo de relajaci´on muy largos y dependencia del protocolo de enfriamiento. Finalmente, en el resto de este trabajo, ⟨(···)⟩denotar´a el promedio termal t´ıpico y (···) denotar´a el promedio sobre el desorden (congelado). 1.1 Modelos de vidrios de esp´ın Se han desarrollado muchos modelos de vidrios de esp´ın para modelizar los sistemas reales, con diferentes formas de afrontar el problema. Vamos a describir brevemente aqu´ı algunos de ellos. El primer tipo de modelo que uno puede estudiar es un sistema que reproduzca el vidrio de esp´ın experimental
CAP´ ITULO 1. INTRODUCCI ´ ON 9 q P(q) qEA Figure 1.3: Representaci´on esquem´atica de la distribuci´on del overlap en la fase paramagn´etica. matriz era de la forma de la Ec. (1.29) ˆ Q0−step = 0q0 ... q00 (1.29) Sin embargo, vimos que esta soluci´on es incorrecta porque, en concreto, la entrop´ıa era negativa, por lo que debemos proponer un nuevo Anzatz. El primer paso consiste en crear n/m grupos de m1r´eplicas cada uno, y tome Qab el valor q1si aybpertenecen al mismo grupo y q0si pertenecen a diferentes grupos. Ahora, la matriz se ha roto en n/m1×n/m1bloques, cada uno de tama˜no m1×m1. Llamemos a esta matriz la matriz de paso 1 (1-step), y en la Ec. (1.30) se muestra un ejemplo de una ˆ Q1−step t´ıpica. ˆ Q1−step = m1 z }| { 0q1 ... q10 q0. . . q0 q0 0q1 ... q10 . . . q0 . . .. . ..... . . q0q0. . . 0q1 ... q10 (1.30)
10 CAP´ ITULO 1. INTRODUCCI ´ ON Con este Ansatz el resultado es mejor que en la soluci´on sim´etrica pero sigue siendo incorrecta (la entrop´ıa sigue siendo negativa pero menor), por lo que se ha de probar un nuevo Ansatz, el segundo paso. Ahora, dividiremos cada grupo en m1/m2×m1/m2bloques, cada uno de tama˜no m2×m2, donde m1ym2son a´un n´umeros enteros. Esta matriz de paso 2 (2-steps) tiene la forma de la Ec. (1.31). m1 z }| { m2 z }| { 0q2 ... q20 . . . q1 . . . ... . . . q1. . . 0q2 ... q20 . . . q0 . . ..... . . q0. . . 0q2 ... q20 . . . q1 . . . ... . . . q1. . . 0q2 ... q20 1111 (1.31) pero la soluci´on sigue siendo incorrecta. Sin embargo, cuanto m´as pasos realizamos mejor es la soluci´on, por lo que si se repite estos pasos infinitas veces, se hallar´a la soluci´on correcta. En este caso los n´umeros enteros mise convierten en una variable continua x∈(0,1) y todos los qmse convierten en una funci´on continua q(x). Por tanto, una matriz de Parisi se pude escribir como Q= (0, q(x)), donde el primer t´ermino es el valor de la diagonal de la matriz (en los ejemplos anteriores era siempre cero porque est´abamos trabajando en ausencia de campo magn´etico externo) y el segundo t´ermino es el valor del resto de elementos de la matriz. En presencia de un campo magn´etico externo el valor de los sitios de la diagonal deja de ser cero, por lo que una matriz de Parisi general es Q= (q, q(x)). Antes de calcular la soluci´on con este Ansatz, vamos a demostrar c´omo trabajar con este tipo de matrices. La traza de la matriz es bastante f´acil de calcular trQ=nq (1.32) Para calcular las siguientes cantidades, las calcularemos primero en un paso finito y despu´es en el l´ımite de infinitos pasos (∞-step). n ∑ a,b Qab =n[q+∑(mi−mi+1)qi]→nq −∫1 n q(x)dx (1.33)
CAP´ ITULO 1. INTRODUCCI ´ ON 11 n ∑ a,b Ql ab =n[ql+∑(mi−mi+1)ql i]→nql−∫1 n ql(x)dx (1.34) Finalmente, calcularemos el producto de dos matrices de Parisi, A= (a, a(x)) yB=(b, b(x)). El resultado es la matriz AB =C= (c, c(x)) donde c=ab −⟨ab⟩(1.35) c(x) = na(x)b(x)+[a−⟨a⟩]b(x) + [b−⟨b⟩]a(x) −∫x n [a(x)−a(y)] [b(x)−b(y)] dy (1.36) y con ⟨a⟩=∫1 n a(x)dx (1.37) Ahora, ya podemos calcular los t´erminos relevantes de la energ´ıa libre, Ec. (1.28), cerca del punto cr´ıtico sin campo magn´etico externo, que se expresa como G(q). G(q) = lim n→0 1 2n[θtrQ2−1 3trQ3−1 6∑ a,b (Qab)4](1.38) El t´ermino cuadr´atico se calcula usando las Ecs. (1.35) y (1.32) trQ2=−n∫1 n q2(x)dx (1.39) El t´ermino cu´artico se calcula usando la Ec. (1.34) ∑ a,b Q4 ab =−n∫1 n q4(x)dx (1.40) Finalmente el t´ermino c´ubico, que es el m´as complicado, se calcula utilizando las Ecs. (1.35), (1.36) y (1.32) trQ3=n[∫1 n xq3(x)dx + 3 ∫1 n dxq(x)∫x n q2(y)dy](1.41) Substituyendo las Ecs. (1.39), (1.40) y (1.41) en la Ec. (1.38) y tomando el l´ımite n→0, se encuentra que la energ´ıa libre G(q) = 1 2∫1 0 dx [|θ|q2(x) + 1 6q4(x)−1 3xq3(x)−q(x)∫x 0 q(y)dy](1.42)
12 CAP´ ITULO 1. INTRODUCCI ´ ON Ahora, la ecuaci´on de punto silla se puede escribir como δG δq(x)= 0 (1.43) y haciendo la derivada funciona, se halla 2|θ|q(x) + 2 3q3(x)−xq(x)−2q∫1 x q(y)dy −∫x 0 xq(x)dx (1.44) y derivando con respecto a xse encuentra |θ|+q2(x)−xq(x)−∫1 x q(y)dy = 0 (1.45) y derivando de nuevo se obtiene finalmente q(x) = x 2or dq dx = 0 (1.46) La soluci´on es q(x) = x/2 para valores peque˜nos de xyq(x) = qmax constante para valores grandes de x(n´otese que si la soluci´on fuera q(x) = q0en todo el rango de x∈(0,1), se recobrar´ıa la soluci´on sim´etrica en las r´eplicas). Sea x1el punto en el que cambia el comportamiento de la soluci´on. Como la soluci´on tiene que ser continua, 2qmax =x1y substituyendo en la Ec. (1.45) se halla que qmax =|θ|+O(θ2) (1.47) En la Figura 1.4 se puede observar esta soluci´on. N´otese que si existe un campo magn´etico externo no nulo, de acuerdo a la Ref. [21], existe otra zona plana para valores peque˜nos de xcon valor qmin(h) = 3 4[h2 J2]2 3 (1.48) Estudiaremos ahora la funci´on de distribuci´on del overlap. En general, se puede escribir que P(q) = 1 n(n−1) ∑ a=b δ(Qab −q) (1.49) Substituyendo Qab con una matriz de Parisi, se halla que P(q) = 1 n(n−1)n[(n−m1)δ(q−q0)+(m1−m2)δ(q−q2) + . . . ] →−1 n−1∫1 n δ[q−q(x)]dx (1.50)
CAP´ ITULO 1. INTRODUCCI ´ ON 13 x q x1 qmax qmin(h) Figure 1.4: Representaci´on esquem´atica de la soluci´on hallada para RSB. La l´ınea de puntos es la zona plana a bajas temperaturas en presencia de un campo magn´etico externo. Finalmente, tomando el l´ımite se llega a P(q) = dx(q) dq (1.51) donde x(q) es la funci´on inversa de q(x). N´otese que en esta soluci´on, P(q) tiene una funci´on delta de Dirac en q=qmax y no es nula en (0, qmax). En la Figura 1.5, se puede observar una representaci´on esquem´atica de este resultado. q P(q) qEA Figure 1.5: Representaci´on esquem´atica de la distribuci´on del overlap de la soluci´on RSB. Finalmente, como corolario, si se calcula la distribuci´on del overlap de
14 CAP´ ITULO 1. INTRODUCCI ´ ON tres r´eplicas, se halla que P(q1, q2, q3) = 1 2P(q1)x(q1)δ(q1−q2)δ(q1−q3) +1 2[P(q1)P(q2)θ(q1−q2)δ(q2−q3) +P(q1)P(q3)θ(q3−q1)δ(q1−q2) +P(q2)P(q3)θ(q2−q3)δ(q3−q1)] (1.52) P(q1, q2, q3) ´unicamente no se anula cuando los tres overlaps son iguales o cuando lo son dos y el tercero es mayor que ambos. Por tanto, los overlaps se organizan siguiendo las normas de un espacio ultram´etrico. 1.3.2 Modelo de droplet La teor´ıa de los droplets fue desarrollada por Bray y Moore [22, 23] usando el grupo de renormalizaci´on de Migdal-Kadanoff [24, 25], y desde un punto de vista fenomenol´ogico por Fisher y Huse [26, 27, 28]. En este caso se trabaja con un Hamiltoniano con interacciones de corto alcance. Un droplet es una regi´on compacta en la que los espines est´an invertidos. La distribuci´on de probabilidad de la energ´ıa libre de un droplet es P[∆F(L)] = 1 Lyf(∆F Ly)(1.53) Calculemos ahora la funci´on de correlaci´on [28] G(rij) = [⟨σiσj⟩−⟨σi⟩⟨σj⟩]2(1.54) AT= 0, esta funci´on de correlaci´on tiende a cero. Sin embargo, a temperatura T≪1 G(rij)∝P[∆F(rij]≃P[0] (1.55) por lo que G(rij)∝1 ryand ξ→ ∞ (1.56) Ahora, si se elige una funci´on de correlaci´on un poco diferente, se puede calcular que G1(rij) = ⟨σiσj⟩2−⟨σi⟩2−⟨σj⟩2∼(q2−q2)∼1 ry→0 (1.57) por lo que la distribuci´on del overlap es bastante simple, como se muestra en la Figura 1.6.
CAP´ ITULO 1. INTRODUCCI ´ ON 15 q P(q) qEA -qEA Figure 1.6: Representaci´on esquem´atica de la distribuci´on del overlap en un droplet. Finalmente, estudiaremos el comportamiento de un droplet en presencia de un campo magn´etico externo, comparando c´omo escalan la energ´ıa de la pared de un droplet y el campo [2]. En primer lugar, la energ´ıa de la pared de un dominio escala como Ly, donde ydebe satisfacer la desigualdad y≤D−1 2(1.58) donde Des la dimensi´on del sistema. Por otra parte, el campo externo escala como LD/2. De acuerdo con la Ec. (1.58), y < D/2 para todas las dimensiones, por lo que el campo crece m´as r´apidamente que la energ´ıa de la pared del dominio, por lo que el orden magn´etico no es estable a largas distancias. 1.3.3 Consecuencias La soluci´on RSB en campo medio (interpret´andolo como un modelo con interacciones de corto alcance en dimensi´on infinita) es la soluci´on exacta por encima y en la dimensi´on cr´ıtica superior, DU= 6, mientras que el modelo droplet es la soluci´on exacta en bajas dimensiones. Sin embargo, no se conoce el comportamiento de un sistema realista en tres dimensiones. Afortunadamente, como hemos comprobado en esta secci´on, el comportamiento esperado en cada escenario es completamente diferente. RSB predice una transici´on de fase a un fase de vidrio de esp´ın en presencia de campo magn´etico mientras que el modelo droplet no. Adem´as, la distribuci´on de probabilidad del overlap es bastante diferente en estos escenarios: en RSB, hemos hallado una funci´on
16 CAP´ ITULO 1. INTRODUCCI ´ ON delta de Dirac en qEA y una probabilidad no nula P(q)>0 para 0 < q < qEA, Figura 1.5; mientras que en el modelo droplet hemos hallado s´olo una funci´on delta de Dirac en qEA, Figura 1.6. Esta diferencia en P(q) ser´a bastante ´util para distinguir si el escenario RSB o el droplet es el correcto. De hecho, hay otro escenario intermedio, TNT (trivial-no trivial) pero nos centraremos en los dos primeros. 1.4 Caracter´ısticas de la fase de vidrio de esp´ın En las secciones anteriores hemos explicado varios modelos de vidrios de esp´ın y la transici´on de fase de algunos de ellos. Ahora, presentaremos algunas caracter´ısticas de los vidrios de esp´ın, especialmente de su fase de vidrio de esp´ın. La frustraci´on es una de las principales causas de los tiempos de relajaci´on largos que exhiben estos sistemas en su fase de vidrio de esp´ın, y como consecuencia, la hip´otesis de ergodicidad deja de cumplirse. Adem´as, una herramienta ´util para detectar transiciones de fase, el par´ametro de orden de los vidrios de esp´ın, ser´a presentado. 1.4.1 Broken ergodicity Si se quiere medir un observable en un experimento, se deber´ıa probar con un tiempo de observaci´on mayor que el tiempo de relajaci´on m´as largo del sistema. De esta forma, el sistema puede explorar todo el espacio de fase, y esta medida es equivalente a un promedio estad´ıstico en equilibrio. Este fen´omeno se conoce como ergodicidad. Sin embargo, en algunos sistemas, estos no ocurren, por ejemplo si el tiempo de relajaci´on diverge en el l´ımite termodin´amico (N→ ∞). Estos sistemas se suelen conocer como no erg´odicos. De acuerdo con la soluci´on RSB de Parisi del modelo SK (Secci´on 1.3.1 y Refs. [17, 18, 19, 20]), ´este es el caso de los vidrios de esp´ın y su fase de vidrio de esp´ın en la que existen varios estados puros (fases). En el l´ımite termodin´amico, estos estados (tambi´en llamados valles porque son los m´ınimos de la energ´ıa libre) tienden a tener la misma energ´ıa libre, pero las barreras entre ellos tienden a infinito. Por tanto, el sistema no puede explorar todo el rango de microestados y el sistema se convierte en no erg´odico. Sin embargo, no est´a claro si este es el comportamiento de un vidrio de esp´ın con interacciones de alcance finito.
CAP´ ITULO 1. INTRODUCCI ´ ON 17 1.4.2 Par´ametro de orden Un par´ametro de orden, si existe, es una herramienta com´un y ´util para estudiar las transiciones de fase. Se trata de un observable cuyo comportamiento en cada fase es diferente. Un ejemplo de par´ametro de orden es el overlap definido en la Secci´on 1.3. Su distribuci´on de probabilidad tiene s´olo una delta de Dirac en q= 0 en la fase paramagn´etica (Figura 1.3), pero en la fase de vidrio de esp´ın tiene dos deltas de Dirac en el escenario droplet (Figura 1.6) o dos deltas de Dirac y una parte continua en el escenario RSB de Parisi (Figura 1.5). Respecto a los vidrios de esp´ın, el primer par´ametro de orden fue propuesto por Edwards y Anderson [13], definido como qEA = lim t→∞ lim N→∞⟨σi(t0)σit0+t)⟩(1.59) donde el promedio t´ermico corre sobre un conjunto de distintos valores de t0. Como comentamos m´as arriba, en el l´ımite termodin´amico, las barreras entre las diferentes fases tienden a infinito, por lo que el sistema no es capaz de cambiarse del valle en el que est´a. Como consecuencia, esta cantidad es una medida de la magnetizaci´on local media, promediada sobre todos los valles. Puede tambi´en escribirse como qEA =∑ a Pa⟨σi⟩2 a(1.60) donde el ´ındice acorre sobre todas las fases y Paes la probabilidad t´ermica Pa=e−βFa ∑ b e−βFb (1.61) Definamos ahora la magnetizaci´on cuadr´atica local media en equilibrio q=⟨σi⟩2(1.62) Calcul´andola s´olo sobre una muestra, se puede escribir como qJ=1 N∑ i⟨σi⟩2=1 N∑ i∑ ab PaPb⟨σi⟩a⟨σi⟩b(1.63) Es tambi´en bastante ´util calcular la correlaci´on entre distintas fases, por lo que definimos el overlap para una sola muestra como qab =1 N∑ i⟨σa i⟩⟨σb i⟩(1.64)
18 CAP´ ITULO 1. INTRODUCCI ´ ON que tiene la propiedad de |qab| ≤ 1. Estudiando la distribuci´on de probabilidad de esta cantidad [29], se halla que P(q) = ⟨δ(q−qab)⟩=∑ ab PaPbδ(q−qab) (1.65) Si solo hay dos fases, como en el modelo droplet,P(q) ser´ıa la suma de dos funciones delta de Dirac (una en −qEA y la otra en qEA), como encontramos antes (Figura 1.6). Por otra parte, si el sistema presenta una rotura de ergodicidad no trivial, como RSB, P(q) tendr´ıa tambi´en una parte continua, como encontramos anteriormente (Figura 1.5). Finalmente, vamos a aplicar el m´etodo de las r´eplicas (Secci´on 1.2) para calcular estas cantidades en el contexto en el que trabajaremos en esta tesis. De acuerdo con Ref. [2], se puede definir qαβ =⟨σα iσβ i⟩(1.66) donde los ´ındices αyβindican un par de r´eplicas distintas (α=β). Por tanto, se puede identificar q= lim n→0 1 n(n−1) ∑ α=β qαβ (1.67) donde nindica el n´umero de r´eplicas. Finalmente, tambi´en se puede identificar el par´ametro de orden de Edwards-Anderson como qEA = max αβ qαβ (1.68) 1.4.3 Frustraci´on La Frustraci´on es una de las mayores contribuciones al tremendamente complicado paisaje de energ´ıa libre que provoca la t´ıpica lenta din´amica de los vidrios de esp´ın. Como ejemplo, vamos a trabajar con un vidrio de esp´ın de tipo Ising de Edwards-Anderson sin campo magn´etico externo. Como se puede observar en la Figura 1.7, como los acoplamientos son aleatorios, ciertas configuraciones de ellos provocan que algunos espines (en la figura el de la esquina inferior derecha) no sean capaces de encontrar la posici´on m´as estable. En nuestro ejemplo, el esp´ın de la esquina inferior derecha tiende a estar apuntando hacia arriba para estar en paralelo con el esp´ın de la esquina inferior izquierda, debido al acoplamiento que hay entre ellos. Sin embargo, tambi´en tiende a estar apuntando hacia abajo para estar antiparalelo al esp´ın de la esquina superior derecha, debido al acoplamiento entre ambos.
CHAPTER 1. INTRODUCTION 25 spin glass model (see Section 1.1). The main quantity is the free energy, that can be computed as FJ=−KBTlog ZJ(1.10) However, one must average over the samples, so F=FJ=−KBTlog ZJ(1.11) The disadvantage of this relation is that averaging a logarithm is quite difficult. The solution is the replica method, based on log Z= lim n→0 Zn−1 n(1.12) Therefore, we have nreplicas of the system and the average over the disorder can be computed as Zn≡Zn J= n ∏ a=1 Z(a) J=∑ {σa i} exp (−β n ∑ a=1 HJ({σa i}))(1.13) In the case of the Edwards-Anderson model this relation becomes in Zn J=∑ {σa i} exp (1 4β2∑ ij J∑ ab σa iσb iσa jσb j)(1.14) so the initial problem of averaging over the disorder have become in a problem of computing ndifferent replicas. Then, one must extend it to a non-integer value of nand take the limit when n→0. From Eq. (1.14) one notices that an effective Hamiltonian, Heff , which depends on the spins of two different replicas can be defined. In fact, every observable that depends on a set of k thermal averages of spins, this observable can be rewritten using kdifferent replicas. 1.3 Approximations of spin glass models with exact solution In Section 1.1, several realistic spin glasses models have been presented. However, the analytic solution of these models is quite difficult, so some approximation should be performed. In this section, we will present two approaches which allow us an analytical computation: the mean field approximation and the droplet model.
26 CHAPTER 1. INTRODUCTION 1.3.1 Sherrington-Kirkpatrick model In 1975, Sherrington and Kirkpatrick [16] propounded a mean field theory based on a model with infinite range interactions. The Hamiltonian of Sherrington-Kirkpatrick (SK) model is H=−1 2∑ i=j Jijσiσj+∑ i hiσi(1.15) where the distribution of the couplings, P(Jij) is Gaussian (the same for every pair of spins) with Jij =J0(1.16) J2 ij =J2 N(1.17) Notice that, comparing this model with the Edwards-Anderson one, Eq. (1.5), SK model is a kind of EA model where every spin interacts with an infinite number of neighbours, thus SK model is usually interpreted as an EA model in infinite dimensions. Symmetric solution (paramagnetic phase) Using the replica method explained in Section 1.2, we will firstly compute the partition function Zn=∑ [σa] exp {−J0β+J2β2 2[n 2(N−1)) −n(n−1) 2] +J0β 2N∑ a(∑ i σa i)2 +J2β2 2N∑ a<b (∑ i σa iσb i)2 (1.18) where the indices aand bruns over the replicas of the system. One needs to avoid quadratic terms, which can be achieved by using the HubbardStratonovich identity. Now, we will change our variables to new ones, Qab and ma, defined as Qab =1 N N ∑ i⟨σa iσb i⟩(1.19) ma=1 N N ∑ i⟨σa i⟩(1.20)
CHAPTER 1. INTRODUCTION 27 As a consequence, an effective partition function can be defined Zeff ≡∑ [σa] exp [(Jβ)2∑ a<b Q2 abσaσb+J0β∑ a maσa](1.21) thus one finally gets Za∝∫[dm] [dQ] exp [−1 2NJ0β∑m2 a−1 2N(Jβ)2∑ a<b Q2 ab +Nlog Zeff] ≡∫[dm] [dQ] exp [−NG(m, Q)] (1.22) where [dm]≡∏ a dma(1.23) [dQ]≡∏ ab dQab (1.24) Eq. (1.22) defines a new function G(m, Q). Let (m0 a, Q0 ab) be the saddle point and let us assume the symmetric solution Ansatz: m0 a≡mand Q0 ab ≡q, that is, all the replicas have the same parameters. Therefore, one can compute the free energy (per spin), f(m, q) which is the function G(m, Q) in the previous Eq. (1.22) f(m, q) = −J2β 4(1−q2)+J0 2m2(1.25) −1 β∫dz √2πe−1 2z2log [2 cosh (Jβ√qz +βh +J0mβ)] and the equilibrium values m=∫dz √2πe−1 2z2tanh (Jβ√qz +βh +J0mβ) (1.26) q=∫dz √2πe−1 2z2tanh2(Jβ√qz +βh +J0mβ) (1.27) It is quite easy to compute that if h= 0, when T > Tfthe only solution is q= 0, but when T < Tf, the observable q(T)= 0 (in fact when T→0, q→1), where Tfis a critical temperature. Thus one has an order parameter. In Figure 1.3, the probability distribution of this order parameter q is represented when T > Tf(the high temperature phase, the paramagnetic one).
28 CHAPTER 1. INTRODUCTION q P(q) qEA Figure 1.3: Schematic representation of the distribution of the overlap in the paramagnetic phase. However, this symmetric solution is not correct, at least at low temperatures, whereas one can assume that it does hold in the paramagnetic phase. The breakdown of the symmetric solution for T < Tfis signaled by a negative value of the entropy at T= 0 and for the appearance of negative eigenvalues in the Hessian matrix. Therefore, one has to compute a new solution at low temperatures that avoids these problems, and this solution will be the Parisi’s Replica Symmetry Breaking (RSB) [17, 18, 19, 20]. Parisi’s Replica Symmetry Breaking Firstly, we will expand the argument of the exponential in Eq. (1.22), so, assuming J0= 0, we can get that G(ˆ Q) = lim n→0 1 n[−1 2τtr(Q2)−1 6tr(Q3)−1 12 ∑ a,b Q4 ab +1 4∑ a=b=c Q2 abQ2 ac −1 8tr( ˆ Q4)]+O(Q5) (1.28) where τ= (Tc−T)/Tcand θ=−τ. One can ignore the two last terms, 1 4∑ a=b=c Q2 abQ2 ac and −1 8tr( ˆ Q4) because they finally get terms O(τ5) or O(τ6) which can be ignored. Once one has defined the free energy in function of the matrix ˆ Q, we will now discuss the Ansatz for this matrix. The first Ansatz one could imagine is the one that we have used in the replica symmetric solution, let us name
CHAPTER 1. INTRODUCTION 29 it the 0-step matrix. Remind that the matrix was like in Eq. (1.29) ˆ Q0−step = 0q0 ... q00 (1.29) However, we saw that this solution is incorrect because, in particular, the entropy was negative, so one can deal with a new Ansatz. The first step consists of creating n/m groups of m1replicas every one, and let Qab be q1if aand bbelong to the same group and q0if they belong to different groups. Now, the matrix has been broken into n/m1×n/m1blocks, every one of size m1×m1. Let us name it the 1-step matrix and in Eq. (1.30), an example of a typical ˆ Q1−step is shown. ˆ Q1−step = m1 z }| { 0q1 ... q10 q0. . . q0 q0 0q1 ... q10 . . . q0 . . .. . ..... . . q0q0. . . 0q1 ... q10 (1.30) With this Ansatz the result is better than in the replica symmetric solution but it is still incorrect (the entropy is also negative but smaller), so a new Ansatz can be tested, the second step. Now one divide every group in m1/m2×m1/m2blocks, every one of size m2×m2, where m1and m2are
30 CHAPTER 1. INTRODUCTION still integers. This 2-steps matrix reads like in Eq. (1.31). m1 z }| { m2 z }| { 0q2 ... q20 . . . q1 . . . ... . . . q1. . . 0q2 ... q20 . . . q0 . . ..... . . q0. . . 0q2 ... q20 . . . q1 . . . ... . . . q1. . . 0q2 ... q20 1111 (1.31) but the solution is also incorrect. However, the more steps one makes, the better the solution is, so if one repeats these steps infinitely times, one will find a correct solution. Then, the integers mitend to a continuous variable x∈(0,1) and all the qmbecome in the continuous function q(x). Therefore, a Parisi’s matrix can be written as Q= (0, q(x)), where the first term is the value in the diagonal of the matrix (in the previous examples it was always zero because we were working in no external magnetic field) and the second term is the value of the rest of the matrix elements. In presence of an external magnetic field, the value of the diagonal sites is not zero, so the general Parisi’s matrix is Q= (q, q(x)). Before computing the solution with this Ansatz, we will show how to work with this kind of matrices. The trace of the matrix is quite easy to compute trQ=nq (1.32) In order to compute the following quantities, we will firstly calculate them in a finite step and later in the limit ∞-step. n ∑ a,b Qab =n[q+∑(mi−mi+1)qi]→nq −∫1 n q(x)dx (1.33) n ∑ a,b Ql ab =n[ql+∑(mi−mi+1)ql i]→nql−∫1 n ql(x)dx (1.34)
CHAPTER 1. INTRODUCTION 31 Finally, we will compute the product of two Parisi’s matrices, A= (a, a(x)) and B=(b, b(x)). The result is the matrix AB =C= (c, c(x)) where c=ab −⟨ab⟩(1.35) c(x) = na(x)b(x)+[a−⟨a⟩]b(x) + [b−⟨b⟩]a(x) −∫x n [a(x)−a(y)] [b(x)−b(y)] dy (1.36) and with ⟨a⟩=∫1 n a(x)dx (1.37) Now, we can compute the relevant terms of the free energy, Eq. (1.28), near the critical point without an external magnetic field, which is denoted as G(q). G(q) = lim n→0 1 2n[θtrQ2−1 3trQ3−1 6∑ a,b (Qab)4](1.38) The quadratic term is computed using Eqs. (1.35) and (1.32) trQ2=−n∫1 n q2(x)dx (1.39) The quartic term is computed using Eq. (1.34) ∑ a,b Q4 ab =−n∫1 n q4(x)dx (1.40) Finally the cubic term, which is the most complicated one, is computed using Eqs. (1.35), (1.36) and (1.32) trQ3=n[∫1 n xq3(x)dx + 3 ∫1 n dxq(x)∫x n q2(y)dy](1.41) Substituting Eqs. (1.39), (1.40) and (1.41) in Eq. (1.38) and evaluating the limit n→0, the free energy reads G(q) = 1 2∫1 0 dx [|θ|q2(x) + 1 6q4(x)−1 3xq3(x)−q(x)∫x 0 q(y)dy](1.42) Now, the saddle point equation can be written as δG δq(x)= 0 (1.43)
32 CHAPTER 1. INTRODUCTION and performing the functional derivative, one obtains 2|θ|q(x) + 2 3q3(x)−xq(x)−2q∫1 x q(y)dy −∫x 0 xq(x)dx (1.44) and differentiating it with respect to xone finds |θ|+q2(x)−xq(x)−∫1 x q(y)dy = 0 (1.45) and differentiating again one finally obtains q(x) = x 2or dq dx = 0 (1.46) The solution is q(x) = x/2 for small xand q(x) = qmax constant for large x (notice that if the solution was q(x) = q0in the whole range of x∈(0,1), the replica symmetric solution would be recovered). Let x1be the point where the change of the behavior of the solutions takes place. As the solution must be continuous, 2qmax =x1and substituting in Eq. (1.45) one finds that qmax =|θ|+O(θ2) (1.47) In Figure 1.4 one can see this solution. Notice that if the external magnetic field does not vanish, according to Ref. [21], there is another plateau at small values xwith value qmin(h) = 3 4[h2 J2]2 3 (1.48) We will now study the overlap distribution function. In general, one can write that P(q) = 1 n(n−1) ∑ a=b δ(Qab −q) (1.49) Substituting Qab with a Parisi’s matrix one finds that P(q) = 1 n(n−1)n[(n−m1)δ(q−q0)+(m1−m2)δ(q−q2) + . . . ] →−1 n−1∫1 n δ[q−q(x)]dx (1.50) Finally, evaluating the limit one gets P(q) = dx(q) dq (1.51)
CHAPTER 1. INTRODUCTION 33 x q x1 qmax qmin(h) Figure 1.4: Schematic representation of the solution found to RSB. Dotted line is the plateau at low temperatures in presence of an external magnetic field. q P(q) qEA Figure 1.5: Schematic representation of the distribution of the overlap the RSB solution. where x(q) is the inverse function of q(x). Notice that in this solution, P(q) has a Dirac’s delta function at q=qmax and it does not vanish in (0, qmax). In Figure 1.5, one can see a schematic representation of this result. Finally, as a corollary, if one computes the distribution of the overlap of
34 CHAPTER 1. INTRODUCTION three replicas one finds that P(q1, q2, q3) = 1 2P(q1)x(q1)δ(q1−q2)δ(q1−q3) +1 2[P(q1)P(q2)θ(q1−q2)δ(q2−q3) +P(q1)P(q3)θ(q3−q1)δ(q1−q2) +P(q2)P(q3)θ(q2−q3)δ(q3−q1)] (1.52) P(q1, q2, q3) does not vanish only when the three overlaps are equal or when two of them are equal and the third one is bigger than them. Hence, the overlaps organize with the rules of an ultrametric space. 1.3.2 Droplet Model The theory of the droplets was developed by Bray and Moore [22, 23] using Migdal-Kadanoff renormalization group [24, 25], and from a phenomenological point of view by Fisher and Huse [26, 27, 28]. In this case one works with a Hamiltonian with short range interactions. A droplet is a compact region of reversed spins. The probability distribution of the free energy of a droplet is P[∆F(L)] = 1 Lyf(∆F Ly)(1.53) We will now compute the correlation function [28] G(rij) = [⟨σiσj⟩−⟨σi⟩⟨σj⟩]2(1.54) At T= 0, this correlation function tends to zero. However at a temperature T≪1 G(rij)∝P[∆F(rij]≃P[0] (1.55) hence, G(rij)∝1 ryand ξ→ ∞ (1.56) Now, if one chooses a bit different correlation function, one can compute that G1(rij) = ⟨σiσj⟩2−⟨σi⟩2−⟨σj⟩2∼(q2−q2)∼1 ry→0 (1.57) so the distribution of the overlap is quite simple, as it is shown in Figure 1.6. Finally, we will study the behaviour of a droplet in presence of an external magnetic field, comparing how the energy of the wall of a droplet and the
Chapter 2 Potts 2.1 Preliminary study The Disordered Potts Glass Model (DPM) has been extremely studied because of its interesting characteristics. In particular, mean field model exhibits a dynamic phase transition for a given number of states p, which makes this model quite useful to study supercooled liquids and glasses. Besides, this model does not have any inversion symmetry (σi→σi) whereas Ising-like models in absence of a magnetic field do, so DPM has been used to modelize systems without inversion symmetry, for example orientational (or quadrupolar) glasses [30, 31], like ortho-hydrogen, and mixed crystals [32, 33], like (KCN)x(KBr)1−x. Finally, DPM plays a similar role as the pure Potts model does in the study of ferromagnets, it allows us to study several kinds of phase transitions, a first and second order thermodynamic phase transition and a dynamic one, just varying the value of the number of states (p) in our simulations. 2.1.1 Mean field analysis In 1985, Gross, Kanter and Sompolinsky [34] studied the mean field theory of the DPM (see also Refs. [2, 35, 36] for a more detailed explanation). The mean field Hamiltonian of DPM is H=−1 2∑ i=j Jijδσiσj(2.1) where pis the number of states that a given Potts spin σican take. The (quenched) couplings, Jij, are Gaussian-distributed random variables with 41
42 CHAPTER 2. POTTS mean J0/N. The order parameter can be defined as qrr′=(⟨δσir⟩− 1 p)(⟨δσir′⟩− 1 p)(2.2) which has the symmetry [2] qrr′=q(δrr′−1 p)(2.3) In the replica method (Section 1.2), qbecomes a matrix and can be expressed as Qαβ =⟨δσασβ⟩− 1 p(2.4) where αand βare replica indices and the thermal average is computed with the replicated Hamiltonian. Now, we can compute the free energy (similarly to what we did in Section 1.3.1) near the critical temperature f(Q) = lim n→0 p−1 2n[θtr (Q2)−1 3tr (Q3)−p−2 6∑ αβ Q3 αβ +y(p) 6∑ αβ Q4 αβ](2.5) where, remind, θ= (T−Tc)/Tc. If we compare this result with the analogous Eq. (1.38) found in Section 1.3.1 to the Sherrington-Kirkpatrick (SK) model, we notice two main differences. In Eq. (2.5) we have two cubic terms instead only one: ∑ αβ Q3 αβ does not vanish because in DPM the symmetry under inversion of the spins does not hold. The other main difference is that in Eq. (2.5) the coefficient of the quartic term is not constant but depends on p. Let p∗the value of pwhere y(p) changes its sign. It is negative for p<p∗ and positive for p > p∗with p∗∼2.8 (see Ref. [34]). In fact y(2) = −1, so if p= 2, the SK model is recovered. Firstly, we will study this model in the region p < p∗where y(p) is negative. If we assume that a continuous Parisi’s solution q(x) holds, the solution would be q(x) = −1 4y(p)[2x−(p−2)] or dq dx = 0 (2.6) Due to the fact that 0 ≤q(x)≤1, the solution presents two plateaus joint by a straight line of slope 1/[2y(p)] (see Figure 2.1a). Therefore, the probability
CHAPTER 2. POTTS 43 distribution of the order parameter qhas two Dirac’s deltas, one at q0= 0 and the another one at q1, and it does not vanish in the (0, q1) region (see Figure 2.1b for more details). Notice that if p= 2, we have the same solution as in the SK model. Since q(x) must also be a non-decreasing function, this solution is only correct if y(p)<0, which agrees with our assumption, and it is incorrect for p > p∗. x q (a) q(x) q P(q) (b) P(q) Figure 2.1: Schematic representation of the solution found for p < p∗. Secondly, we will study the region where p > p∗. We will use the one step Replica Symmetry Breaking (RSB) Ansatz, Eq. (1.30), which is enough to solve the system. Let nbe the total number of replicas and m1the number of replicas of every of the n/m1groups. The expression of the Eq. (1.34) in the one step Ansatz becomes (we assume that the terms of the diagonal of the matrix ˆ Qare 0) ∑ α,β Ql α,β =n[(m1−1) ql 1+ (n−m1)ql 0](2.7) Computing tr (Q3) is a bit more tricky, but after some algebra, it can be expressed as tr (Q3)=n{(m1−1) (m1−2) q3 1+(n m1−1)[3m1(m1−1) q2 0q1 +(n m1−2)m2 1q3 0]} (2.8) Substituting Eqs. (2.7) and (2.8) in Eq. (2.5) and evaluating the limit, one
44 CHAPTER 2. POTTS gets f(q) = p−1 2{θ[(m1−1) q2 1−m1q2 0]−1 3[(m1−1) (m1−2) q3 1 −3m1(m1−1) q2 0q1+ 2m2 1q3 0]−p−2 6[(m1−1) q3 1−m1q3 0] +y(p) 6[(m1−1) q4 1−m1q4 0]}(2.9) The saddle point equations can be expressed as 0 = ∂f ∂q0 =p−1 2{−2θm1q0+ 2m1(m1−1) q0q1−2m2 1q2 0(2.10) +p−2 2m1q2 0−2y(p) 3m1q3 0} 0 = ∂f ∂q1 =p−1 2{2θ(m1−1) q1−(m1−1) (m1−2) q2 1(2.11) +m1(m1−1) q2 0−p−2 2(m1−1) q2 1+2y(p) 3(m1−1) q3 1} 0 = ∂f ∂m1 =p−1 2{θ(q2 1−q2 0)−1 3[(m1−1) (m1−2) q3 1−3m1q2 0q1(2.12) −3m1(m1−1) q2 0q1+ 2m1q3 0]−p−2 6(q3 1−q3 0)+y(p) 6(q4 1−q4 0)} Taken into account that q0≤q1and neglecting the quartic term, the solution of these equations, q, is a step function q={0 if x < x0 2θ p−4if x > x0(2.13) where x0is the parameter m x0≡m1=p−2 2(2.14) In Figure 2.2a one can see a schematic representation of this function. The probability distribution of the order parameter qcan be observed in Figure 2.2b. Nevertheless, this solution also becomes incorrect in the region p > 4, where a discontinuous transitions appears. Therefore, the approximation
CHAPTER 2. POTTS 45 x q x0 (a) q(x) q P(q) (b) P(q) Figure 2.2: Schematic representation of the solution found for p > p∗and T2< T < Tc. used to compute Eq. (2.5) is not valid and the previous demonstration does not hold. However, Eq. (2.5) can still be used in the limit ϵ≡p−4→0 because the discontinuity is small enough. In this situation one finds that at the critical temperature, Tc, the value of the order parameter above the discontinuous jump is q(1) ∝p−4 and the position of the jump as temperature tends to the critical one is x0(T→T− c)→1. Whereas, if one solves the full problem when p→ ∞, the value of qabove the discontinuous jump is q(1) = 1 and the position of the jump is x0=T/Tc. Finally, Gross, Kanter and Sompolinsky [34] also demonstrated that the system undergoes a second phase transition at a temperature T2< Tc, because the previous solution has a negative entropy at T= 0 for every finite p > p∗. Using an expansion of the free energy up to fifth order terms in ˆ Q, they demonstrated that this phase transition is a continuous one, so q(x) has a continuous part in a range of x, as can be observed in Figure 2.3a. In Figure 2.3b the probability distribution of the order parameter is plotted. 2.1.2 Glass phase transition The glass phase transition was firstly studied in the framework of the supercooled liquids. This area deal with amorphous solids like the glass of windows. Many reviews of this topic have been written, but we will focused on Refs. [37, 38]. If one cools fast enough a liquid, it would not become solid at its melting temperature, Tm, and it would remain liquid even at temperatures below that temperature Tm: this is a supercooled liquid. However, the lower the temperature, the slower the dynamic of the system, that is, the relaxation time exhibits an extremely growth (several orders of magnitude) in a short
46 CHAPTER 2. POTTS x q (a) q(x) q P(q) (b) P(q) Figure 2.3: Schematic representation of the solution found for p > p∗and T < T2< Tc. range of temperature. In fact, at a temperature low enough the relaxation time is so long that the system is not able to explore the whole phase space in the time that a typical experiment spend, so the system becomes non-ergodic. This behavior defines a kind of dynamical phase transition, where the phase at low temperature is the so-called glass phase. To compute the temperature at which the phase transition takes place, Tg, one needs to establish a criterion to determine whether an experimental time is long enough to characterize a glass phase. This maximum experimental time is usually fixed at 102−103s. With this definition, the viscosity where the glass phase transition happens can be computed [38]: η(Tg) = 1013 Poise.(2.15) Different liquids have not the same evolution of the viscosity. Some of them, strong liquids, have a fast evolution, linear versus Tg/T, for example SiO2. Other liquids, fragile liquids, have a far slower evolution at high temperatures, for example o-terphenyl. In Figure 2.4, this behaviour can be observed. Notice that, where Tg/T = 1, all liquids have the same evolution due to the definition of Tg, Eq. (2.15). This definition of the phase transition and Tgseems to be a mathematical trick without any physical meaning. In fact Tgdepends (weakly) on the cooling protocol of the experiment. However, this is not the case thanks to some characteristics of these systems, like the so-called two steps relaxation. Let C(t1, t2) be a general defined two times correlation function C(t1, t2) = 1 N∑ i⟨ϕi(t1)ϕi(t2)⟩(2.16) where ϕis an observable that depends on the particle (in liquids) or on the spin (in spin glasses) which stays in the position i. An example of this kind
CHAPTER 2. POTTS 47 Figure 2.4: Evolution of the viscosity of several liquids. Notice that, although the evolution is different, all of them reach the same value of the viscosity. This figure is the famous Angell plot, from Ref. [45] of two times correlation function in spin glasses is the one defined in Eq. (6.2), where the observable is the spin itself. At equilibrium, the two times correlation function, Eq. (2.16), does just depend on the difference of these times t=t2−t1, so Eq. (2.16) can be rewritten as C(t) = 1 N∑ i⟨ϕi(t)ϕi(0)⟩(2.17) At high temperature, C(t) decrease as an exponential function C(t) = Aexp(−t/τ) (2.18) However, this behaviour does not hold at temperatures near Tg, where a plateau in the relaxation of C(t) appears, that is, C(t) decreases and reaches a first plateau and later it resumes the decreasing. The length of this plateau depends on the temperature and appears continuously as Tdecreases, so this phase transition is usually called a continuous transition. Nevertheless, if one focuses on the value of C(t) on the plateau, one has a discontinuous behaviour as Tdecreases. Some characteristics of these supercooled liquids and their glass transition seems to be quite similar to properties of spin glasses in their spin
48 CHAPTER 2. POTTS glass phase, such as the extremely long relaxation time. Besides, spin glasses without reversal spins symmetry, like DPM or p-spin model [39], undergo a discontinuous phase transition (in fact, several phase transitions actually happens, some of them continuous), as it is shown in the previous section (Section 2.1.1) for DPM. However, this phase transition seems to be a first order one, at least in mean field analysis. Kirkpatrick, Wolynes and Thirumalai [40, 41, 42, 43, 44] performed an in-depth study of this relation between spin glasses (they specially worked with Potts glass model) and the glass transition of the supercooled liquids. For example, they found [41] that the correlation function exhibit a plateau, a behaviour similar to the two steps relaxation. 2.1.3 Previous results Brangian, Kob and Binder [46, 47, 48] performed a complete study of the tenstate infinite range DPM and Gaussian couplings with a negative mean. They checked whether this model presents the dynamical and static phase transitions that mean field theory predicts in the thermodynamic limit. Therefore they performed simulations of several system sizes, up to N= 2560 spins. They simulated 500 samples for the smallest system and between 20 and 50 for the largest one. They found strong finite size effects, although their simulations suggest the existence of both static and dynamical transition. Therefore a finite system behaves, at least qualitatively, similarly as in the thermodynamic limit. Regarding the more realistic short range models, Brangian, Kob and Binder [49, 50] also studied the three dimensional ten-state short range model, although they focused on a bit different model from the one we will study in this chapter H=−∑ ⟨i,j⟩ Jij (pδσiσj−1)(2.19) where Gaussian and bimodal couplings were studied, both with a negative mean (J0<0). Systems sizes up to L= 16 have been simulated, with up to 100 samples and 108Monte Carlo steps (MCS). For both probability distributions of the quenched couplings, they did not find any sign of the existence neither the static nor the dynamical transition predicted in mean field theory, so the behavior of the short range systems would be extremely different to the infinite range ones. Lee, Katzgraber and Young [55] also studied the short range DPM (in fact they studied the same model we will study in this chapter). They performed simulations of the four dimensional three-state DPM and the three
CHAPTER 2. POTTS 49 dimensional threeand ten-state DPM, all of them with Gaussian couplings. Besides, the three dimensional three-state model was also studied with bimodal couplings. In all of the three-state models, the probability distribution of the couplings was chosen with a vanishing mean, J0= 0, but in the tenstate model the mean was chosen negative, J0=−1. They simulated in the three dimensional three-state systems of size up to L= 12 performing ∼107MCS and 352 samples in the Gaussian probability distribution and 550 samples in the bimodal one (in both cases, more samples were simulated in smaller lattice sizes). They found a clear phase transition with both probability distributions. However, in the three dimensional ten-state Gaussian DPM1, they did not find any sign of phase transition, which supports the previous result of Brangian, Kob and Binder [49, 50]. Finally, the Janus Collaboration [63] studied the three dimensional four states DPM with binary quenched couplings (with vanishing mean J0= 0). They used a prototype board of Janus (see Section 7.6, Appendix A and Refs. [224, 225, 226]) to perform their simulations. The simulations performed were far longer than in previous works. The statistic achieved was astonishing: the largest lattice size simulated was L= 16 with 1000 samples and 8 ×109MCS every one. Therefore, they were more confident that the system was completely thermalized. They found a clear phase transition to a spin glass phase. Moreover, no sign of a ferromagnetic phase transition was found. Therefore, an in-depth study of the three dimensional DPM with p > 4 states is quite interesting. In this work, we have used the full Janus machine which allows us to perform long simulations with far more statistic than previous works, so we are far more confident that our simulations are completely thermalized. Therefore, it seems that the results obtained in this study are more reliable, although some of them do not agree with previous works, which performed shorter simulations. However we did not manage to thermalize systems of L= 16 lattice size for p≥5 nor even relevant lattice sizes for p= 8. There are several open questions in this model that we will try to understand better with this work. For example whether a phase transition to a ferromagnetic transition exists or characterize the spin glass transition, its order and the behavior of βcwith p. 1The largest lattice size simulated was L= 12 with 343 samples and ∼104MCS (in smaller sizes, 1000 samples were simulated with up to ∼105MCS).
50 CHAPTER 2. POTTS Critical Behavior of Three-Dimensional Disordered Potts Models with Many States R. Alvarez Ba˜nos, A. Cruz, L. A. Fernandez, A. Gordillo-Guerrero, J. M. Gil-Narvion, M. Guidetti, A. Maiorano, F. Mantovani, E. Marinari, V. Martin-Mayor, J. Monforte-Garcia, A. Mu˜noz Sudupe, D. Navarro, G. Parisi, S. Perez-Gaviro, J. J. Ruiz-Lorenzo, B. Seoane, S. F. Schifano, A. Tarancon, R. Tripiccione and D. Yllanes. Published in J. Stat. Mech. P05002 (2010). 2.2 Introduction The three dimensional (3D) disordered Potts model (DPM) is an important system, that could help in clarifying a number of open and crucial questions. The first issue that comes to the mind is the possibility of understanding the glass transition, since this is a very challenging problem. On more general grounds, it is very interesting to try and qualify the behavior of the system when the number of states pbecomes large: here we should see the paradigm of a “hard”, first order like transition but, as we will discuss in the following, only sometimes this turns out to be clear (see for example the set of large scale, very accurate numerical simulations of Ref. [51], dealing with a model slightly different from the one defined here). In such a difficult situation extensive numerical simulations are more than welcome, and the Janus supercomputer [52, 53], optimized for studying spin glasses, reaches its peak performances when analyzing lattice regular systems based on variables that can take a finite, small number of values: disordered Potts models fit very well these requirements. Using the computational power of Janus we have been able to consistently thermalize the DPM with p= 5 and 6 on 3D(simple cubic) lattice systems with periodic boundary conditions and size up to L= 12. Bringing these systems to thermal equilibrium becomes increasingly harder with increasing number of states: it has been impossible for us, even by using a large amount of time of Janus (that for these problems performs, as we discuss better in the following, as thousands of PC processors), to get a significant, unbiased number of samples thermalized, and reliable measurements of physical quantities, for p≥5 on a L= 16 lattice. Our results lead us to the claim that the critical behavior of the DPM with a large number of states pis very subtle, and if pis larger than, say, 5, numerical simulations could easily give misleading hints. The numerical results that we will discuss in the following lead us to believe that the spin glass transition gets stronger with increasing number of states p: a theoretical
CHAPTER 2. POTTS 57 1 2 3 4 5 1021031041051061071081091010 ξ Monte Carlo Sweep L = 4 , β = 9.5 L = 6 , β = 9.5 L = 8 , β = 6.5 L = 12 , β = 5.5 Figure 2.5: Log-binning thermalization test for p= 5. For all data points the point size is bigger than the corresponding error bar. later times. Sample to sample fluctuations of τint are very large: in figure 2.8 we plot τint for all our samples with p= 5, L= 8. In order to be on the safe side we have increased the number of MCS, by continuing the numerical simulation for a further extent, in all samples where our estimate of τint was bigger than the length of the simulation divided by a constant c(c= 20 for L= 8 and c= 15 for L= 12, where achieving thermalization is much more difficult). 5 5In the p= 5, L= 8 case for 2442 samples we have run a simulation of total extent η= 4 ×108MCSs, while for 5 samples η= 8 ×108MCSs, and for 1 sample η= 1.6×109 MCS. In the p= 5, L= 12 case for 2382 samples η= 6 ×109MCSs, for 54 samples η= 1.2×1010, for 8 samples η= 2.4×1010, and for 7 samples η= 4.8×1010 MCS. In the p= 6, L= 8 case: for 1263 samples η= 109MCSs, for 8 samples η= 2 ×109and for 9 samples η= 4 ×109. In the p= 6, L= 12 case for 1173 samples η= 6 ×1010 MCSs, for 17 samples η= 1.2×1011 MCSs and for 6 samples η= 2.4×1011 MCSs.
58 CHAPTER 2. POTTS 1 2 3 4 5 1021031041051061071081091010 1011 ξ Monte Carlo Sweep L = 4 , β = 9.8 L = 6 , β = 9.65 L = 8 , β = 7.5 L = 12 , β = 6.5 Figure 2.6: As in figure 2.5, but p= 6. 2.5.2 Critical temperature and critical exponents Our analysis of the critical exponents of the system has been based on the quotient method [56, 62]: by using the averaged value of a given observable Omeasured in lattices of different sizes, we can estimate its leading critical exponent xO, ⟨O(β)⟩ ≈ |β−βc|−xO.(2.33) By considering two systems on lattices of linear sizes Land sL respectively one has that [56, 62] ⟨O(β, sL)⟩ ⟨O(β, L)⟩=sxO/ν +O(L−ω),(2.34) where νis the critical exponent of the correlation length and ωis the exponent of the leading-order scaling-corrections [56]. We use the operators ∂βξ, from (2.26), and χq, from (2.25) in equation (2.34) to obtain respectively the critical exponents 1 + 1/ν and 2 −ηq. The exponent 2 −ηmis obtained applying eq. (2.34) to the magnetic susceptibility χm, from (2.28).
CHAPTER 2. POTTS 59 Figure 2.7: The autocorrelation function (2.31) for one generic sample (p= 6, L= 8).
60 CHAPTER 2. POTTS 0 2 4 6 8 10 12 14 16 0 500 1000 1500 2000 τint Sample Figure 2.8: Integrated autocorrelation time, τint, for all p= 5, L= 8 samples. τint is in units of blocks of ten measurements, i.e. of 20103MCS. Samples above the green line have been “extended” (see the text for a discussion of this issue).
CHAPTER 2. POTTS 61 To use the quotient method we start estimating the finite-size transition temperature: we do this by looking at the crossing points of the correlation length in lattice units (ξ/L) for various lattice sizes. We have used a cubic spline interpolating procedure to compute both the crossings of ξ/L and its β-derivative (we have followed the approach described in detail in Ref. [63]). 0.1 0.2 0.3 0.4 2 3 4 5 6 7 ξ / L β L = 4 L = 6 L = 8 L =12 0.35 0.4 4.5 5 5.5 Figure 2.9: Overlap correlation length in lattice size units as a function of the inverse temperature βfor L= 4, 6, 8 and 12. Here p= 5. We show in figures 2.9 and 2.10 the behavior of ξ/L as a function of β. The different curves are for different lattice sizes. The crossing points are rather clear in both cases, giving a strong hint of the occurrence of a second order phase transition. At least for p= 5 scaling corrections play a visible role, and the crossing points undergo a small but clear drift towards lower temperatures for increasing lattice sizes. We summarize in tables 2.3 and 2.4 the βvalues of the crossing points for two different pairs of lattice sizes, together with the estimated values of the critical exponents νand ηqthat we obtain using relation (2.34). Since we can only get reliable results on small and medium size lattice we cannot control in full scaling corrections, and a systematic extrapolation to the infinite volume limit is impossible. It is clear however that the effective critical exponents summarized in tables 2.3 and 2.4 do not suggest that
62 CHAPTER 2. POTTS 0.1 0.2 0.3 0.4 2 3 4 5 6 7 ξ / L β L = 4 L = 6 L = 8 L =12 0.35 0.4 6 6.5 Figure 2.10: As in figure 2.9, but p= 6. asymptotically for large volume the system will not be critical (in this case, for example, ηqshould be asymptotically equal to 2): our numerical data clearly support the existence of a finite temperature phase transition. (L1, L2)βcross(L1, L2)ν(L1, L2)ηq(L1, L2)ηm(L1, L2) (4,8) 4.83(5) 0.82(3) 0.13(2) 1.72(2) (6,12) 5.01(4) 0.81(2) 0.16(2) 1.94(2) Table 2.3: Numerical values of our estimates for the crossing point of the curves ξ/L. We give βcross, the thermal critical exponent ν, the anomalous dimension of the overlap ηq, and the anomalous dimension of the magnetization ηm. We take as our best estimates for the critical exponents the one obtained from the lattices with sizes L= 6 and L= 12. For p= 5 βc= 5.01(4) , ν = 0.81(2) , ηq= 0.16(2) ,(2.35)
CHAPTER 2. POTTS 63 (L1, L2)βcross(L1, L2)ν(L1, L2)ηq(L1, L2)ηm(L1, L2) (4,8) 6.30(9) 0.80(2) 0.10(2) 1.453(19) (6,12) 6.26(7) 0.80(4) 0.16(2) 1.971(19) Table 2.4: As in table 2.3, but p= 6. while for p= 6. βc= 6.26(7) , ν = 0.80(4) , ηq= 0.16(2) .(2.36) It is interesting to compare these values with those of other Potts models with a different number of states. In particular we are interested in the value of the critical exponents as a function of the number of states, since we want to characterize the critical behavior of the various models and attempt a prediction of the model’s behavior when the number of states is large. In our particular model and with the (low) values of the temperature that are interesting for us (since we need to get below the critical point) even with the large computational power available to us thanks to Janus the simulation for p= 8, say, on a L= 12 lattice, would require an unavailable amount of CPU time. What is found in the very interesting work of Refs. [51] and [55] is different, since there one is able to thermalize a p= 10 model on a large lattice, and no transition is observed. The model analyzed in these two references [51, 55] is indeed slightly (or maybe, it will turn out, not so slightly) different from the present one, since there Jis negative. It is not clear to us if this difference could explain a quite dramatic discrepancy of the observed behavior, or if, for example, a different (very low) temperature regime should be analyzed to observe relevant phenomena: this is surely an interesting question to clarify, and the fact that the coupling have a negative expectation value, reducing in this way frustration, could turn out to make a difference. 2.5.3 Absence of ferromagnetic ordering in the critical region Our DPM is in principle allowed to undergo a ferromagnetic phase transition (since no symmetry protects it), and at low temperatures could present a spontaneous magnetization, as discussed in Ref. [[63]]. Because of that we have carefully studied the magnetic behavior of the model at low temperatures. We have analyzed both the magnetization and the magnetic susceptibility below the spin glass critical point.
64 CHAPTER 2. POTTS 0 1 2 3 4 5 6 1 2 3 4 5 6 7 8 9 10 χm β L = 4 L = 6 L = 8 L = 12 0 0.1 0.2 0.3 2 4 6 8 10 〈|m|〉 β Figure 2.11: Magnetic susceptibility as a function of βfor L= 4, 6, 8 and 12. Here p= 5. In the paramagnetic phase the magnetization is random in sign, and its absolute value is expected to be proportional to 1/√V. In Figs. 2.11 and 2.12 we check whether ⟨|m|⟩ around the spin glass critical region tends to an asymptotic value for larger lattice size, or not. From the figures we see ⟨|m|⟩ goes to zero in the critical region. Also, we studied the magnetic susceptibility χm=V⟨|m|2⟩which is independent of size. Again in Figs. 2.11 and 2.12 we check that, and we see a non-divergent behavior. This behavior is extremely different from a ferromagnetic phase in which χmdiverges as the volume. Besides, as reported in Sec. 2.5.2 the exponent ηmis close to 2, so we could say that a ferromagnetic-paramagnetic phase transition does not happen in the range of temperatures that we have studied. 2.6 Evolution of critical exponents with p In table 2.5 we summarize the values of the inverse critical temperature and of the thermal and overlap critical exponents for DPM from p= 2 (the Ising,
CHAPTER 2. POTTS 65 0 1 2 3 4 5 6 7 8 9 10 2 3 4 5 6 7 8 9 10 χm β L = 4 L = 6 L = 8 L = 12 0 0.1 0.2 0.3 2 4 6 8 10 〈|m|〉 β Figure 2.12: As in figure 2.11, but p= 6. Edwards-Anderson spin glass) up to p= 6. We also plot these data items in figure 2.13. From table 2.5 and figure 2.13 some results emerge very clearly. First, the inverse critical temperature roughly follows a linear behavior in p, with a slope very close to one. We have added in table 2.5 the ratio (R) between the numerical determinations (in 3d) of βc(p) and their values in the Mean Field (MF) approximation. One can see that the large deviations from the MF prediction occur for large values of p(notice that R > 1 since MF suppresses fluctuations). 6 6In the MF approximation was obtained, using the Hamiltonian [57, 67], H ≡ −p 2∑ i=j Jij δsi,sj, that Tc/J = 1 for p≤4 and (Tc/J)2= 1 + (p−4)2/42 + O((p−4)4) for p > 4. In addition for very large p,Tc/J ≃1 2(p/ log p)1/2. Taking into account the extra pfactor in the Hamiltonian used in the Mean Field and the fact that J=√2d(J2 ij =J2/N, being Nthe number of spins in the MF computation) since we are working in finite dimension (d), we obtain the finite dimension version of the critical βusing the Mean Field approximation: βc=p/√2dfor p≤4 and βc=p √2d(1−(p−4)2/84 + O((p−4)4))
66 CHAPTER 2. POTTS Second, νdecreases monotonically and ηqgrows monotonically with the number of states p. To discuss this behavior it is useful to keep in mind that when using finite size scaling to study a disordered first order phase transition one expects to find [64] ν= 2/D and 2 −ηq=D/2, i.e., in our D= 3 case, ν= 2/3 and ηq= 1/2. These are “effective” exponents, that are a bound to the ones allowed for second order phase transitions. Both sets of values for νand ηqare indeed completely compatible with tending, as pincreases, to those limit values that characterize a first order phase transition. If this turns out, as our numerical data make very plausible to be true, two different scenarios open. The first possibility is that the pstates DPM undergoes a disordered first order phase transition for large enough values of p(just as in the ordered Potts model, that for p≥3 undergoes a first order phase transition), while the second possibility is that the DPM will show a standard second order phase transition for all finite values of p. This is the typical issue that is very difficult to settle with numerical work: an analytical solution of the model with infinite number of states would be very useful as a starting point in order to discriminate between these two possible scenarios. p βcν ηqR 2 (Ref.[[65]]) 1.786(6) 2.39(5) 7−0.366(16)82.187(8) 2 (Ref.[[66]]) 1.804(16) 2.45(15) −0.375(10) 2.209(20) 3 (Ref.[[55]]) 2.653(35) 0.91(2) 0.02(2) 2.17(3) 4 (Ref.[[63]]) 4.000(48) 0.96(8) 0.12(6) 2.45(3) 5 (this paper) 5.010(40) 0.81(2) 0.16(2) 2.51(2) 6 (this paper) 6.262(71) 0.80(4) 0.16(2) 2.69(3) Table 2.5: Critical parameters as a function of p. All data are for binary couplings, with zero expectation value. By Rwe denote the ratio between the critical βin three dimensions and that computed in Mean Field. 2.7 Conclusions In this note we have characterized the critical behavior of the 3DDPM with p= 5 and p= 6, i.e. with a reasonably large number of states. Our numerical for p > 4 (notice the minus signum of the (p−4)2correction); in addition, for large p, one obtains βc≃√2 d(plog p)1/2. Note that in our case √2d≃2.45.
CHAPTER 3. SAMPLE TO SAMPLE FLUCTUATIONS 73 3.1.2 Overlap equivalence The overlap equivalence is the property of a system when every generalized overlap that one can define, using an arbitrary observable, O, qO=1 N∑ i Oi(σa)Oj(σb) (3.23) depends on the usual overlap q. That is, although they both do fluctuate when N→ ∞,qOrestricted to pairs of replicas with a given qdo not fluctuate. Therefore, the usual overlap qcontains all the useful information and, thus, is the complete order parameter3. Separability is a similar property, but using equilibrium configurations instead of real replicas. These two properties are equivalent4. Let Mab be matrices that belong to the set of all matrices computed from the matrix Q(for example ∑ c QacQcb). In 2000, Parisi and Ricci-Tersenghi [74] demonstrated that the overlap equivalence (or separability) implies that ∑ b Qk abMab =∑ b∫dqδ(q−Qab)Qk abMab =∫dqP(q)M(q)qk(3.24) where M(q) is the value that the matrix Mab takes when Qab =qand P(q) = ∑ b δ(q−Qab) (3.25) They setted that ∑ b Qk abMab =∫dqP(q)M(q)qk(3.26) ∑ b Qk abM′ ab =∫dqP(q)M′(q)qk(3.27) ∑ b Qk abMabM′ ab =∫dqP(q)M(q)M′(q)qk(3.28) (3.29) Notice that the probability of the last equation can be written in function of the probabilities of the other two relations P(q)M(q)M′(q) = [P(q)M(q)] [P(q)M′(q)] P(q)(3.30) 3For example, in mean field, the energy overlap satisfies qe=q2. 4See Ref. [74] for more details of this property. Moreover, a detailed analysis of the overlap equivalence is performed in this reference.
74 CHAPTER 3. SAMPLE TO SAMPLE FLUCTUATIONS Choosing two matrices like Mab =∑ c Qk1 acQk2 cb (3.31) M′ ab =∑ c Qk3 acQk4 cb (3.32) and considering all possible values of k, they found that the joint probability P(5) ≡P12,13,32,24,41 can be computed as 3P(5)(q, q1, q2, q3, q4) = δ(q1−q4)δ(q2−q3)P(3)(q, q1, q2) + 2P(3)(q, q1, q2)P(3)(q, q3, q4) P(q)(3.33) where P(3) is defined as P(3) ≡P12,23,31 (3.34) Now, integrating Eq. (3.33) over q, the joint probability P(4) ≡P13,32,24,41 is computed 3P(4)(q1, q2, q3, q4) = 1 2δ(q1−q4)δ(q2−q3) [P(q1)P(q2) + δ(q1−q2)P(q2)] + 2 ∫dqP(3)(q, q1, q2)P(3)(q, q3, q4) P(q)(3.35) This relation is quite useful because P(4)(q1, q2, q3, q4) is, by construction, invariant under permutations of the overlaps, but the right hand term of Eq. (3.35) is not for a generic P(3) function. Therefore this equation enforces hard constrains in P(3). Impose equations like P(4)(qi, qj, qj, qi)−P(4)(qi, qi, qj, qj) = 0 (3.36) P(4)(qi, qi, qj, ql)−P(4)(qi, qj, ql, qi) = 0 (3.37) One can compute admissible P(3). 3.1.3 Replica equivalence and overlap equivalence imply ultrametricity Once we have explained the concepts of replica equivalence and overlap equivalence, we will now demonstrate that if this two concepts hold, then ultrametricity also hold.
CHAPTER 3. SAMPLE TO SAMPLE FLUCTUATIONS 75 In order to compute P(q) and P(3), we assume that overlap can take only a few values, k(we suppose that our results also hold in the continuous case), so these probabilities are just a sum of delta functions P(q) = k ∑ i=1 piδ(q−qi) (3.38) P(3)(qi, qj, ql) = ∑ i,j,l pijlδ(q−qi)δ(q−qj)δ(q−ql) (3.39) where piand pijl are weights, the last one is invariant under permutations of the indices. Obviously all the weights belong to the interval [0,1]. Besides, relations like Eqs. (3.17) and (3.20) (from stochastic stability or replica equivalence) generate new relations between these pijl weights. Remind that in Section 1.3.1, where the existence of ultrametricity is shown in a infinite range model (but using distances defined from the overlaps), Eq. (1.52) tells us that only equilateral and isosceles triangles are allowed (in fact some isosceles triangles are also forbidden, those that do not satisfy this relation, that is, assuming q1< q2< q3, all weights piij vanish if i>j). Therefore, if the scalene terms, pijl with i=j=land the forbidden isosceles terms vanish, ultrametricity holds. In fact, after using the symmetry under permutations of the indices and relations from replica equivalence, only weights from scalene and forbidden isosceles are free parameters, the rest of the parameters can be expressed as a function of them and weights pi. For a given k, there are (k 3)scalene weights and (k 2)forbidden isosceles parameters. Now, we will study, as an example, the case when k= 5. Using Eqs. (3.35) and (3.36) with the overlaps q4and q5, we get 0 = 1 4p5p4+p2 541 p1 +p2 542 p2 +p2 543 p3 +p2 544 p4 +p2 554 p5 −(p551p441 p1 +p552p442 p2 +p553p443 p3 +p554p444 p4 +p555p544 p5)(3.40) Using replica equivalence relations, the allowed isosceles parameter and the equilateral parameters of the previous relation can be written as a function of the forbidden isosceles and scalene parameters p555 =1 2p5(1 + p5)−p554 −p553 −p552 −p551 (3.41) p444 =1 2p4(1 + p4−p5)−p441 −p442 −p443 +p541 +p542 +p543 +p554 (3.42) p544 =1 2p4p5−p541 −p542 −p543 −p554 (3.43)
76 CHAPTER 3. SAMPLE TO SAMPLE FLUCTUATIONS Substituting Eqs. (3.41), (3.42) and (3.43) in Eq. (3.40) we find 0 = 2p1p2p3p4p5E5,4 0+ [2p1p2p3p5p554 +p1p2p3p4(3p5−2p555)] (p543 +p542 +p541)+4p1p2p3p5(p543p542 +p543p541 +p542p541) + [2p2p3p5(p1+p4)] p2 541 + [2p1p3p5(p2+p4)] p2 542 + [2p1p2p5(p3+p4)] p2 543 (3.44) where in E5,4 0we include all the terms independent of the scalene parameters. It can be written as E5,4 0=p1p4p5 4[2p551 p1p5(1−2p554 p4p5)(1−2p441 p1p4)+(1−2p551 p1p5)4p554p441 p1p2 4p5] +p2p4p5 4[2p552 p2p5(1−2p554 p4p5)(1−2p442 p2p4)+(1−2p552 p2p5)4p554p442 p2p2 4p5] +p3p4p5 4[2p553 p3p5(1−2p554 p4p5)(1−2p443 p3p4)+(1−2p553 p3p5)4p554p443 p32p2 4p5](3.45) The terms with the form 2piij pipj (3.46) belong to [0,1] due to the fact that all the weights are positive. Therefore, it is obvious that E4,3 0is non-negative. Thus, all of the terms of Eq. (3.44) are also non-negative, so in order to satisfy the equation, all the scalene parameters must vanish. Repeating this method with other pairs of overlaps, one finds that all the scalene parameters do vanish p543 =p542 =p541 =p532 =p531 =p521 =p432 =p431 =p421 =p321 = 0 (3.47) The following step is to study the isosceles parameters (in this case in k= 4), although they are a bit more tricky. Using Eqs. (3.35) and (3.37) with the overlaps q1,q2and q3(assuming q1< q2< q3< q4), we get 0 = p1p2p3 4[2p331 p1p3(1−2p332 p2p3)(1−2p221 p1p2)+(1−2p331 p1p3)2p332 p2p3 2p221 p1p2](3.48) Repeating this method one finds other three similar relations. All of them imply that three of the forbidden isosceles parameters vanish p331 =p441 =p442 = 0 (3.49) and one of the rest p332,p221 or p443 do also vanish. Therefore, two of the forbidden isosceles parameters do not vanish and ultrametricity is violated. Fortunately, for a general kthere are (k 2)∼k2forbidden isosceles parameters and k/2 of these parameters violate ultrametricity. As kgrows proportion of isosceles parameters which violate ultrametricity decreases and in the limit k→tends to 0. We can conclude that if one assumes that replica equivalence and overlap equivalence hold, all the scalene and forbidden isosceles parameters vanish and, thus, ultrametricity also holds.
CHAPTER 3. SAMPLE TO SAMPLE FLUCTUATIONS 77 3.1.4 Ultrametricity in short range models Finally, in 1996 ´ I˜niguez, Parisi and Ruiz-Lorenzo [75] demonstrated that if ultrametricity holds in a short range spin glass model, Eq. (1.52) of ultrametricity in mean field is recovered. Let Hbe this short range spin glass model H=−∑ ⟨ij⟩ Jijσiσj(3.50) which is invariant under permutations of replicas. Then, the general expression of the joint probability P12,13,23 is P12,13,23(q12, q13, q23) = A(q12)δ(q12 −q13)δ(q12 −q23) +B(q12, q13)θ(q12 −q13)δ(q13 −q23) +B(q13, q23)θ(q13 −q23)δ(q23 −q12) +B(q23, q12)θ(q23 −q12)δ(q12 −q13) (3.51) Moreover, the two replicas probability P12,13 can be computed from Eq. (3.51) integrating over q23, so P12,13(q12, q13) = [A(q12) + ∫∞ q12 dq23B(q13, q23)]δ(q12 −q13) +B(q12, q13) (3.52) Integrating again, now over q13, the one replica probability distribution is computed P(q12) = A(q12) + ∫∞ q12 dq13B(q12, q13) + ∫∞ −∞ dq13B(q12.q13) (3.53) Using Eq. (3.17) and a little algebra, one finds that A(q12) = ∫q12 −∞ dq13B(q12, q13) (3.54) B(q12, q13)=2(∫∞ −∞ dq23B(q12, q23))(∫∞ −∞ dq23B(q13, q23))(3.55) Besides, using Eqs. (3.53) and (3.54) one gets P(q12) = 2 ∫∞ −∞ dq23B(q12, q23) (3.56) Finally, taking into account Eqs. (3.56), (3.54) and (3.55) one finds A(q12) = 1 2x(q12)P(q12) (3.57) B(q12, q13) = 1 2P(q12)P(q13) (3.58)
78 CHAPTER 3. SAMPLE TO SAMPLE FLUCTUATIONS Substituting them in Eq. (3.51), Eq. (1.52) is recovered. So, if ultrametricity holds in finite dimensional spin glasses, it will be the same kind of ultrametricity as obtained in mean field
CHAPTER 3. SAMPLE TO SAMPLE FLUCTUATIONS 79 Sample-to-sample fluctuations of the overlap distributions in the three-dimensional Edwards-Anderson spin glass R. A. Ba˜nos, A. Cruz, L. A. Fernandez, J. M. Gil-Narvion, A. Gordillo-Guerrero, M. Guidetti, D. I˜niguez, A. Maiorano, F. Mantovani, E. Marinari, V. Martin-Mayor, J. Monforte-Garcia, A. Mu˜noz-Sudupe, D. Navarro, G. Parisi, S. Perez-Gaviro, F. Ricci-Tersenghi, J. J. Ruiz-Lorenzo, S. F. Schifano, B. Seoane, A. Taranc´on, R. Tripiccione and D. Yllanes. Published in Phys. Rev. B 84, 174209 (2011). 3.2 Introduction Spin glasses are model glassy systems which have been studied for decades and have become a paradigm for a broad class of scientific applications. They not only provide a mathematical model for disordered alloys and their striking low-temperature properties (slow dynamics, age-dependent response), but they have also been the test-ground for new ideas in the study of other complex systems, such as structural glasses, colloids, econophysics, and combinatorial optimization models. The non-trivial phase-space structure of the mean-field solution to spin glasses [76, 77, 78] encodes many properties of glassy behavior. Whether the predictions of the mean-field solutions correctly describe the properties of finite-range spin-glass models (and of their experimental counterpart materials) is a long-debated question. The Droplet Model describes the spin glass phase in terms of a unique state (apart from a global inversion symmetry) and predicts a (super-universal) coarsening dynamics for the off-equilibrium regime. [79] Moreover, there is no spin glass transition in presence of any external magnetic field. On the other side, the Replica Symmetry Breaking scenario [78, 80], based on the mean field prediction, describes a complex free-energy landscape and a non-trivial order parameter distribution in the thermodynamic limit; the dynamics is critical at all temperatures in the spin-glass phase. The spin glass transition temperature is finite also in presence of small magnetic fields; the search for the de Almeida-Thouless line Tc(h) is the purpose of many numerical experiments (see, for example, Ref. [81]). From the theoretical perspective, the last decade has seen a strong advance in the understanding of the properties of the mean-field solution: its correctness has been rigorously proved thanks to the introduction of new concepts and tools, like stochastic stability or replica and overlap equivalence [82, 83, 84, 85, 86]. Besides, numerical simulation has been the methodology
80 CHAPTER 3. SAMPLE TO SAMPLE FLUCTUATIONS of choice when investigating finite-range spin glasses, even if the computational approach is severely plagued by the intrinsic properties (slow convergence to equilibrium, slowly growing correlation lengths) of the simulated system’s (thermo)dynamics. In this respect, a Moore-law-sustained improvement in performance of devices for numerical computation and new emerging technologies in the last years has allowed for very fast-running implementation of standard simulation techniques. By means of the non-conventional computer Janus [87] we have been able to collect high-quality statistics of equilibrium configurations of three-dimensional Edwards-Anderson spin glasses, well beyond what would have been possible on conventional PC clusters. Theoretical predictions and Janus numerical data have been compared in detail in Refs. [88] and [89]. One of the main results presented therein is that equilibrium properties at a given finite length scale correspond to out-of-equilibrium properties at a given finite time scale. On experimentally accessible scales (order 104seconds waiting times corresponding to order 102lattice sizes) the Replica Symmetry Breaking picture turns out as the only relevant effective theory. Theories in which some of the fundamental ingredients of the mean-field solutions are lacking (overlap equivalence in the TNT model [90], ultrametricity in the Droplet Model) show inconsistencies when their predictions are compared to the observed behavior. In this work we reconsider the analysis of the huge amount of data at our disposal, focusing on the sample-to-sample fluctuations of the distribution of the overlap order parameter. The assumptions of the mean-field theory allow us to make predictions on the joint probabilities of overlaps among many real replicas which can be tested against numerical data for the three-dimensional Edwards-Anderson model. The structure of the paper is as follows: in section 3.3 we give some details on the considered spin-glass model and the performed Monte Carlo simulations. In the subsequent section we first recall some fundamental concepts such as stochastic stability, ultrametricity, replica and overlap equivalence and some predictions on the joint overlap probability densities, and then present a detailed comparison with numerical data. In section 3.5 we show how finite-size numerical overlap distributions compare to the mean-field prediction in which finite-size effects are appropriately introduced. We finally present our conclusions in the last section.
CHAPTER 3. SAMPLE TO SAMPLE FLUCTUATIONS 81 3.3 Monte Carlo Simulations 3.3.1 The Model We consider the Edwards-Anderson model [91] in three dimensions, with Ising spin variables σi=±1 and binary random quenched couplings Jij =±1. Each spin, set on the nodes of a cubic lattice of size V=L3(Lbeing the lattice size), interacts only with its nearest neighbors under periodic boundary conditions. The Hamiltonian is: H=−∑ ⟨i,j⟩ Jijσiσj,(3.59) where the sum extends over nearest-neighbor lattice sites. In what follows we are dealing mainly with measures of the spin overlap qab =1 L3∑ i σa iσb i,(3.60) where aand bare replica indices, and the sample-dependent frequencies NJ(qab) with which we estimate the overlap probability distribution PJ(q) of each sample (we indicate one-sample quantities by the subscript J): PJ(qab) = ⟨δ(qab −1 L3∑ i σa iσb i)⟩,(3.61) where ⟨(···)⟩is a thermal average. In what follows (···) denotes average over disorder. 3.3.2 Numerical Simulations We present an analysis of overlap probability distributions computed on equilibrium configurations of the three-dimensional Edwards-Anderson model defined in Eq. (3.59). We computed the configurations by means of an intensive Monte Carlo simulation on the Janus supercomputer. Full details of these simulation can be found in Ref. [89].For easy reference, we summarize the parameters of our simulations in Table 3.1. In order to reach such low temperature values, it has been crucial to tailor the simulation time, on a sample-by-sample basis, through a careful study of the temperature randomwalk dynamics along the parallel tempering simulation.
82 CHAPTER 3. SAMPLE TO SAMPLE FLUCTUATIONS L Tmin Tmax NTNS 8 0.150 1.575 10 4000 16 0.479 1.575 16 4000 24 0.625 1.600 28 4000 32 0.703 1.549 34 1000 Table 3.1: A summary of parameters of the simulations we have used in this work. For each lattice size, L, we considered NSsamples, with four independent real replicas per sample. For the Parallel Tempering algorithm, NTtemperatures were used between Tmin and Tmax, uniformly distributed in that range (except in the case of L= 8, in which we have 7 temperatures uniformly distributed between 0.435 and 1.575 plus the 3 temperatures 0.150, 0.245 and 0.340). Our MCS consisted of 10 Heat-Bath sweeps followed by 1 Parallel Tempering update. More detailed information regarding these simulations can be found in Ref. [89]. 3.4 Replica equivalence and ultrametricity The Sherrington-Kirkpatrick (SK) model [76] is the mean-field counterpart of model (3.59). It is defined by the Hamiltonian H=∑ i=j Jijσiσj,(3.62) where the sum now extends to all pairs of NIsing spins and the couplings Jij are independent and identically-distributed random variables extracted from a Gaussian or a bimodal distribution with variance 1/N. The quenched average of the thermodynamic potential may be performed by rewriting the n-replicated partition function in terms of an n×noverlap matrix Qa,b for which the saddle-point approximation gives the self-consistency equation Qab =⟨σaσb⟩,(3.63) where the average ⟨(···)⟩involves an effective single-site Hamiltonian in which Qa,b couples the replicas. The thermodynamics of model (3.62) is recovered in the limit n→0, after averaging over all possible permutations of replicas. The overlap probability distribution P(q) is defined in terms of such an averaging procedure: for any function of the overlap f(q), one has that ∫dqa,bP(qa,b)f(qa,b) = lim n→0 1 n!∑ p f(Qp(a),p(b)),(3.64)
CHAPTER 3. SAMPLE TO SAMPLE FLUCTUATIONS 89 such effects cannot be clearly told by comparing only the smallest lattices considered, L= 16 and L= 24. At T= 0.57Tc, size-dependent effects are strong even for L= 16,24 (see Fig. 3.5, bottom). Having data from four independent replicas per sample, we have access to the joint probability of two independent overlaps. According to Eq. (3.66) the quantity P(q12, q34) P(q34)−2 3P(q12) = P(q12|q34)−2 3P(q12),(3.85) (where P(·|·) denotes conditional probability) when plotted versus q12, should be a delta function in q34. This quantity is shown for L= 32, T∼0.64Tcand two values of q34 in the top plot of Fig. 3.6 and reveals a clear peak around q34. At high q12 values there is a small excess in the probability P(q12)P(q34), so the difference in Eq. (3.85) becomes negative. As one sees in Fig. 3.6 this happens at values q12 ≳qEA, i.e., in a region of atypically large overlaps that should vanish in the thermodynamical limit. The size dependence for the quantity in Eq. (3.85) is not easy to quantify from the data: as one can see in Fig. 3.6 (bottom) for a particular choice of q34, the peak height tends to increase with L(at least for T∼0.75Tc), but in a very slow way, making extrapolations in the L→ ∞ limit practically impossible. Despite this, we note that the negative peaks get narrower as the system size increases: we expect then that this effect will disappear at larger system sizes. We conclude this section commenting the asymptotic behavior of the cumulative probability ΠC q(z), Eq. (3.81). The small-zdecay is clearly a power law (see top plot in Fig. 3.7), but the best fit exponent is significantly different from the estimate obtained by integrating the overlap distribution P(q). Fig. 3.7 shows a comparison of the exponent x(q) obtained by the two methods, for some lattice sizes, many cut-off values qand two temperatures, T∼0.64Tcand T∼0.57Tc. Although the differences seem to decrease by increasing the lattice size, the trend is very slow and even not in a clear direction for some values of the cutoff q. Again, the only conclusion that can be drawn is that the finite-size effects are large, even for L= 32, and safe extrapolations in the L→ ∞ limit cannot be done. A closer inspection of the data reported in Fig. 3.7 reveals that the difference between the two data sets is roughly a constant, and this difference becomes extremely important in the limit of small q, where one would expect both measurements of x(q) to approach zero. Contrary to expectations, the x(q) estimated from the data of ΠC qseems to remain non-zero even in the q→0 limit. A possible explanation for this observation comes from the fact that the delta peaks in the PJ(q) get broader for systems of finite size. Indeed, in the thermodynamic limit, one would expect PJ(q) to be the sum of
90 CHAPTER 3. SAMPLE TO SAMPLE FLUCTUATIONS delta functions centered on overlap values extracted from the average distribution P∞(q): if this expectation is true, then the value for XJ(q) is nothing but the probability of having a peak at an overlap value smaller than qand this is exactly x(q). However, if the delta peaks acquire a non-zero width ∆ due to finite-size effects, then for q < ∆ the overlap probability distribution close to the origin PJ(0) may be affected by broad peaks centered on overlaps larger than q, which should not count in the thermodynamical limit. If this explanation is correct, then the limit q→0 for the data shown in Fig. 3.7 (bottom) obtained from ΠC qshould give a rough estimate, in the large L limit, for the peak width ∆ (see data in Table 3.3 and discussion below). 3.5 The order parameter distribution We now compare the P(q) obtained in numerical simulations of the threedimensional Edwards-Anderson model (3.59) to the prediction obtained by smoothly introducing controlled finite-size effects on a mean-field-like distribution consisting in a delta function centered in q=qEA and a continuous tail down to q= 0 (a similar analysis has been carried out for long-range spin-glass models, see Ref. [107]). On the positive qaxis one has P∞(q) = e P(q)Θ(qEA −q) + [1 −x∞(qEA)]δ(q−qEA),(3.86) x∞(qEA) = ∫qEA 0 dq e P(q).(3.87) It is convenient to introduce the effective field htrough q= tanh (h) (3.88) and consider its distribution P∞(h) = P∞(q(h))dq(h) dh =dq(h) dh e P(q(h))Θ(hEA −h) + [1 −x∞(qEA)]δ(h−hEA),(3.89) x∞(qEA) = ∫hEA 0 dh e P(h),(3.90) being clear that qEA = tanh (hEA). This change of variable smooths the constraint on the fluctuations of qnear the extremes of the distribution.
CHAPTER 3. SAMPLE TO SAMPLE FLUCTUATIONS 91 L T/TcqEA x∞(qEA) ∆ 32 0.75 0.663(19) 0.91(13) 0.0923(80) 0.64 0.7319(30) 0.828(28) 0.1015(30) 24 0.75 0.69674(72) 1.0000(3) 0.10618(84) 0.64 0.7625(27) 0.876(24) 0.1182(24) 0.57 0.7954(24) 0.842(25) 0.1216(32) 16 0.75 0.73780(73) 1.000031(7) 0.1443(10) 0.64 0.809(16) 1.00(14) 0.150(11) 0.57 0.8210(41) 0.811(49) 0.1683(51) 8 0.75 0.8250(21) 1.000001(9) 0.2872(37) 0.57 0.886(18) 0.95(18) 0.296(28) L T/Tcα γ χ2/d.o.f. 32 0.75 1.92(34) 11.2(1.2) 20/97 0.64 0.93(44) 7.7(1.0) 38/103 24 0.75 2.04(21) 9.68(55) 45/101 0.64 0.95(21) 6.88(41) 69/107 0.57 0.75(17) 5.62(30) 88/110 16 0.75 1.76(16) 5.14(31) 77/107 0.64 0.45(21) 4.50(52) 133/113 0.57 0.53(19) 3.37(40) 161/115 8 0.75 0.73(22) 2.02(34) 501/121 0.57 0.49(16) 1.36(17) 466/123 Table 3.3: Results of the fitting procedure of Eq. (3.94) on numerical P(q) data, with kernel exponent k= 2.5 (see Eq. (3.91)). All errors on parameters are jackknife estimates. We used the symbol χ2in the table to denote the sum of squares of residuals, which is not a true chi-square estimator as the values of P(q) at different qare mutually correlated. In a finite-size system the thermodynamical distribution P∞(h) will be modified, mainly by the fact that delta functions become distributions with non-zero widths. Remember that, in the thermodynamical limit, we expect the distribution PJ(h) for any given sample to be the sum of delta functions. A simple way to take into account the spreading of the delta functions due to finite-size effects is to introduce a symmetric convolution kernel G(k) ∆(h−h′)≡Cexp [−(|h−h′|/∆)k],(3.91) where Cis a normalizing constant and the spreading parameter ∆ is assumed
92 CHAPTER 3. SAMPLE TO SAMPLE FLUCTUATIONS not to depend on h,5while it should have a clear dependence on the system size, such that limL→∞ ∆ = 0. The parameter k, to be varied in the interval [2,3], is introduced in order to consider convolutions different from the Gaussian case (k= 2). In order to obtain an analytic expression for the finite size distribution PL(h)≡∫dh′P∞(h′) + P∞(−h′) 2G(k) ∆(h−h′),(3.92) we assume the following form for the continuous part of the distribution e P(h)≡e P(q(h))dq(h) dh =e P(0)(1 + αh2+γh4),(3.93) where e P(0) = e P(0) = P∞(0), αand γare free parameters to be inferred from the data. The final result is PL(h) = [1 −x∞(qEA)]G(k) ∆(h−hEA) + G(k) ∆(h+hEA) 2 +e P(0) ∫hEA −hEA dz [1 + αz2+γz4]G(k) ∆(h−z) (3.94) where x∞(qEA) = 2 e P(0)[hEA +αh2 EA/3 + γh5 EA/5]. We let α,γ,qEA and ∆ vary in a fitting procedure to P(q) Monte Carlo data; values of e P(0) are fixed to the Monte Carlo values PMC (0). The choice of the exponent kin the convolution kernel is crucial. We varied kin the interval [2,3]. The Gaussian convolution k= 2 turned out to be the worst choice in such interval, giving rise to unphysical negative weights for the delta function contributions, i.e., 1 −x∞(qEA)<0. We obtained very good results with the choice k= 2.5. Fit parameters are reported in Table 3.3 for some lattice sizes and temperatures, while Fig. 3.8 shows comparison between Monte Carlo P(q) and the relative fitting curve. Although the fitting curves interpolate nicely the numerical P(q), some of the fitting parameters may look strange: in particular qEA is a bit larger than the peak location and x∞(qEA)≃1 (for example, in the L= 32 data the difference is around 2%). It is worth remembering that in the solution of the SK model at low temperatures the continuous part P(q) has a divergence for q→q− EA, which can easily dominate the delta function in finite-size systems (where delta peaks are broadened). Indeed, by increasing the system size, qEA seems 5This introduces a q-dependent spread, as the Jacobian of the transformation (3.88) stretches the distribution at high qvalues.
CHAPTER 3. SAMPLE TO SAMPLE FLUCTUATIONS 93 to move towards the location of the peak maximum and x∞(qEA) becomes smaller than 1. In order to make a stronger test of the above fitting procedure, we have used the fit parameters in Table 3.3 to derive the finite-size conditional probability PL(q|q′) = PL(q, q′)/PL(q′) (3.95) applying the convolution kernel G(k) ∆(h−h′) to the L=∞joint probability given by the Ghirlanda-Guerra relation, r.h.s of Eq.(3.66). Fig. 3.8 shows a comparison between our extrapolated PL(q12|q34 =q0) and the Monte Carlo data for L= 32, T= 0.64Tcand three values of q0: the agreement is very good at any value of q0, especially considering that the fitting parameters were previously fixed by interpolating the unconditional overlap distribution PL(q). 3.6 Conclusions We performed a direct inspection of stochastic stability and ultrametricity properties on the sample-to-sample fluctuations of the overlap probability densities obtained by large-scale Monte Carlo simulations of the threedimensional Edwards-Anderson model. We found small but still sizeable deviations from the prediction of the Ghirlanda-Guerra relations but a clear tendency towards improvement of agreement with increasing system size. Large fluctuations make it difficult to draw any definitive conclusion on the analysis of the ultrametric relation (3.78) when taking into account data for the largest lattice size. In addition, critical effects show up at T∼0.75Tc. Considering that for a stochastically stable system overlap equivalence is enough to infer ultrametricity, the results presented here support and integrate the analyses and claims of Refs [88], [89] and [100], in which the authors reported strong evidence of overlap equivalence. We also turned our attention to the shape of the overlap probability distribution, showing that finite-size PL(q) and PL(q, q′) compare well with mean-field (infinite-size) predictions, modified by finite-size effects that only make delta functions broad.
94 CHAPTER 3. SAMPLE TO SAMPLE FLUCTUATIONS 0.01 0.1 1 0.4 0.6 10-6 10-4 10-2 100 0.01 0.1 1 XT q L=24, T ≈ 0.64Tc ABC AAB AAA 0.01 0.1 1 0.4 0.6 10-6 10-4 10-2 100 0.01 0.1 1 XT q L=32, T ≈ 0.64Tc ABC AAB AAA Figure 3.1: (Color online) The quantity XTas defined in the text, as a function of qfor lattice size L= 24 (top) and L= 32 (bottom) at temperature T≃0.64Tc. Insets show a magnified view of the region q∼0.6 (log-log plot). Plots show data for XTcomputed only with triplets of independent configurations (ABC), with triplets in which two configurations belong to the same Monte Carlo history (AAB), and triplets in which all configurations come from the same Monte Carlo history (AAA). No significant difference shows up as long as we take enough uncorrelated configurations from the same replica.
CHAPTER 3. SAMPLE TO SAMPLE FLUCTUATIONS 95 0 0.25 0.5 0.75 1 0 0.2 0.4 0.6 0.8 1 X2 (X1+2X1 2)/3 L=32 L=24 L=16 f(x)=x 0 0.25 0.5 0.75 1 0 0.2 0.4 0.6 0.8 1 X2/X1 X1 L=32 L=24 L=16 f(x)=(1+2.*x)/3 0 3 10-4 6 10-4 9 10-4 0 0.2 0.4 0.6 0.8 1 K2 X1 L=32 L=24 L=16 Figure 3.2: (Color online) Top: X2as a function of the corresponding polynomial in X1(Eq. (3.76)). The straight line is the theoretical prediction (unit slope). Center: the ratio X2/X1as a function of X1, where the straight line is the theoretical prediction. Bottom: the squared difference K2= [X2−(X1+ 2X2 1)/3]2as function of X1. Data refer to T∼0.64Tc
96 CHAPTER 3. SAMPLE TO SAMPLE FLUCTUATIONS 0 0.25 0.5 0.75 1 0 0.2 0.4 0.6 0.8 1 X3 (2XT+2X1+6X1 2+5X1 3)/15 L=32 L=24 L=16 0 10-3 2 10-3 0 0.2 0.4 0.6 0.8 1 K3 X1 L=32 L=24 L=16 Figure 3.3: (Color online) Data at T∼0.64Tc. Top: X3as a function of the corresponding polynomial in X1and XT(Eq. (3.77)). The straight line is the theoretical prediction (unit slope). Bottom: the squared difference K3= [X3−(2XT+ 2X1+ 6X2 1+ 5X3 1)/15]2as function of X1,T= 0.64Tc. Lines connecting points are only a guide to the eye.
CHAPTER 3. SAMPLE TO SAMPLE FLUCTUATIONS 97 0 0.8 10-3 1.6 10-3 2.4 10-3 0 0.2 0.4 0.6 0.8 1 [XT-X1 2]2 X1 L=32 L=24 L=16 0 0.8 10-3 1.6 10-3 2.4 10-3 0 0.2 0.4 0.6 0.8 1 K3 u X1 L=32 L=24 L=16 Figure 3.4: (Color online) Top: The squared difference [XT−X2 1]2as a function of X1. Bottom: the quantity Ku 3= [X3−(2X1+ 8X2 1+ 5X3 1)/15]2 as a function of X1. All data for T∼0.64Tcand for lattice sizes L= 16,24,32. The lines connecting the data points are only intended as a guide to the eye.
98 CHAPTER 3. SAMPLE TO SAMPLE FLUCTUATIONS 0 10-3 2 10-3 3 10-3 0.1 0.5 0.9 0.1 0.5 0.9 0 10-3 2 10-3 3 10-3 [XT-X1 2]2K3 u X1 T≈0.75Tc L=32 L=24 L=16 L=8 0 10-3 2 10-3 3 10-3 0.1 0.5 0.9 0.1 0.5 0.9 0 10-3 2 10-3 3 10-3 [XT-X1 2]2K3 u X1 T≈0.57Tc L=24 L=16 L=8 Figure 3.5: (Color online) Square difference [XT−X2 1]2(left) and the quantity Ku 3= [X3−(2X1+ 8X2 1+ 5X3 1)/15]2(right) as a function of X1. Top: for T= 0.75Tcand L= 8,16,24,32. Bottom: for T= 0.57Tcand L= 8,16,24.
CHAPTER 4. MICROSCOPIC DYNAMICS OF THE 3D SPIN GLASS IN PRESENCE OF A MAGNETIC FIELD 105 q P(q) qEA Figure 4.3: Probability distribution of qin presence of an external magnetic field in Droplet scenario. 4.1.1 Experimental results Regarding experiments in spin glasses, we will focus on Fe0.5Mn0.5TiO3which is supposed to be a short range Ising spin glass1. This material was studied by J¨onsson et al [108] in a large range of external magnetic fields, up to h= 20000 Oe and no phase transition was reported. They studied the decay of the overlap to control whether the phase transition happened. In figure 4.4, their results about this observable are shown. From the classic Ogielski’s paper [113], it is known that a decay like q(t)∼1/txindicates the onset of a spin glass phase. For h= 1000 Oe in the Figure 4.4, one can observe that the behavior of q(t) is almost a power law, and for h= 300 Oe this behavior is quite clear. Therefore this property suggests us that for a smaller external magnetic field, the spin glass phase transition might have been detected. Moreover, we will show the data of this experiment with those of a one dimensional long range model, KAC model [109, 110] with ρ= 1.5 [110] in Figure 4.5. This model roughly corresponds to the four dimensional short range model. In this figure, one can observe that the critical field in four dimensions is near h= 1000 Oe and critical field decrease with the dimensionality, which supports the previous deduction that experimentalist should try smaller magnetic fields (h≤1000 Oe) to detect a spin glass phase transition in real samples (D= 3). Besides, an AT line was found in Heisenberg spin glasses [111] (remind that Droplet scenario states that this line should not exist even in Heisenberg 1Fe0.5Mn0.5TiO3is an Ising like spin glass whereas AgMn at 2.5% and CdCr1.7IN0.3S4 are Heisenberg like spin glasses.
106 CHAPTER 4. MICROSCOPIC DYNAMICS OF THE 3D SPIN GLASS IN PRESENCE OF A MAGNETIC FIELD Figure 4.4: Behavior for the dynamical overlap, q(τ) (which is proportional to the quantity plotted in the y-axis), as a function of time for different magnetic fields. Figure from Ref. [108]. Figure 4.5: Relative decrease of Tc(h)/Tc(0) with increase field for ρ= 1.5 and h= 0, 0.1, 0.15 and 0.2 versus the relative decrease of χ∗(ZFC susceptibility). Figure from Ref. [110]. Experimental data from Fe0.5Mn0.5TiO3, see Ref. [108]. spin glasses). The same authors also studied Ising-like samples (FeNiPBAl) [112] and reached the same conclusion that we stated here, the magnetic fields used in experiments are too high to see a spin glass phase transition in Ising spin glasses. To sum up, one can conclude that experimental data suggest us that the spin glass phase transition may take place for h < 1000 Oe for the Ising Universality class.
CHAPTER 4. MICROSCOPIC DYNAMICS OF THE 3D SPIN GLASS IN PRESENCE OF A MAGNETIC FIELD 107 4.1.2 Analytical approaches From the theoretical point of view, one can study the replicated Hamiltonian but only above the critical temperature, that is, in the paramagnetic phase. This Hamiltonian becomes [114], in terms of the original overlap field Qab (remind that Qaa = 0), H=1 4∑(∇Qab)2+1 4r∑Q2 ab −1 6w∑QabQbcQca −1 8u∑QabQbcQcdQda +1 4x∑Q2 abQ2 ac −1 8y∑Q4 ab −1 2h2∑Qab +O(Q5, h2Q2).(4.1) where w,u,x, and yare positive couplings. Let Qbe the minima and qab the fluctuations around that minima, then one can write that Qab =Q+qab. In presence of a magnetic field, h, the minima must satisfy [114] rQ + 2wQ2−3uQ3+ 2xQ3−yQ3=h2(4.2) Taking the limit n→0 and neglecting higher orders of qab, the Hamiltonian, Eq. (4.1), becomes in H=1 4∑(∇qab)2+1 4(r+uQ2+ 2xQ2−3yQ2)∑q2 ab −1 2Q(w−uQ −2xQ)∑ a=b qabqac −1 4uQ2∑ a=b=c=d qabqcd −1 6w∑qabqbcqca −1 2uQ2∑ a=c qabqbcqcd +xQ ∑qabq2 ac −1 2yQ ∑q3 ab +O(q4) (4.3) Therefore, the starting field theory is a ϕ3theory with an upper critical dimension DU= 6.Hence, the external magnetic field does not change DU. The critical exponents in six dimensions are ν= 1/2, β= 1 and η= 0, so one has the same critical behavior as in the h= 0 case at D≥6. Notice that we can recover the usual ϕ3field theory for an Ising spin glass in absence of a magnetic field by putting Q= 0: H=1 4∑(∇qab)2+1 4r∑q2 ab −1 6w∑qabqbcqca .(4.4)
108 CHAPTER 4. MICROSCOPIC DYNAMICS OF THE 3D SPIN GLASS IN PRESENCE OF A MAGNETIC FIELD In the case h= 0 (that is, Q= 0) only one propagator exists. However, in presence of an external magnetic field one has three different types of propagators: longitudinal (L), anomalous (A), and replicon (R). Diagonalizing the quadratic term in Eq. (4.3) for finite n, one obtains [115] GL=G1+ 2(n−2)G2+1 2(n−2)(n−2)G3=1 p2+r−2wQ(n−2) ,(4.5) GA=G1+ (n−4)G2−(n−3)3G3=1 p2+r−wQ(n−4) ,(4.6) GR=G1−2G2+G3=1 p2+r+ 2wQ ,(4.7) being nthe number of replicas and pthe momentum. In terms of the original spin variables, G1,G2and G3can be written G1(x) = ⟨sisi+x⟩2,(4.8) G2(x) = ⟨sisi+x⟩⟨si⟩,(4.9) G3(x) = ⟨si⟩2⟨si+x⟩2.(4.10) Notice that if one sets n= 0, Eqs. (4.5) and (4.6) are identical: GA(p) = GL(p) = 1 p2+r+ 4wQ .(4.11) Therefore, one actually has the replicon mode and two degenerated modes (anomalous and longitudinal). Of course, if one sets Q= 0, then the standard propagator is recovered: GA(p) = GL(p) = GR(p) = 1 p2+r.(4.12) In the standard mean field picture, the de Almeida-Thouless line is defined by imposing that only the replicon mode is massless, that is G−1 R(p= 0) = 0, but the other two degenerated modes are massive. In other words, G−1 L(p= 0) = G−1 A(p= 0) >0. Bray and Roberts [114] projected the original theory, Eq. (4.3), into the replicon subspace, using the behavior of the propagators in Mean Field and the degeneration of L and A modes. This is equivalent to setting the longitudinal and anomalous masses to infinity (one can write
CHAPTER 4. MICROSCOPIC DYNAMICS OF THE 3D SPIN GLASS IN PRESENCE OF A MAGNETIC FIELD 109 m2 L∝G−1 L(p= 0) and analogously for the other two modes). The final projected Hamiltonian is: H=1 4∑(∇Rab)2+1 4˜r∑R2 ab −1 6w1∑RabRbcRca −1 6w2∑R3 ab (4.13) They studied this projected Hamiltonian using a perturbative renormalization group and, at the order of the perturbation theory they used, no fixed points were found. Therefore a new strategy has been developed: •One needs to avoid the degeneration between the anomalous and longitudinal propagators (or masses). So we will work with non zero nand at the very end of the computation, nwill be set to 0. •Due to the fact that the degeneration between the modes L and A has been broken, one can try to explore more exotic scenarios for the Almeida-Thouless line, like mR=mA= 0 and mL>0. •The starting Hamiltonian should be the most general cubic Hamiltonian compatible with symmetry, extending the interacting cubic Hamiltonian from four couplings (as in Bray and Roberts [114]). This strategy has been devised and followed by De Dominicis and Temesvari in reference [116]. Their Hamiltonian (H=H(2) +H(3)) reads H(2) =1 2∑ p[(1 2p2+m1)∑ αβ ϕαβ pϕαβ −p+m2∑ αβγ ϕαγ pϕβγ −p+m3∑ αβγδ ϕαβ pϕγδ −p] (4.14) H(3) =−1 6√N∑′ p1p2p3[w1∑ αβγ ϕαβ p1ϕβγ p2ϕγα p3+w2∑ αβ ϕαβ p1ϕαβ p2ϕαβ p3(4.15) +w3∑ αβγ ϕαβ p1ϕαβ p2ϕαγ p3+w4∑ αβγδ ϕαβ p1ϕαβ p2ϕγδ p3+w5∑ αβγδ ϕαβ p1ϕαγ p2ϕβδ p3 +w6∑ αβγδ ϕαβ p1ϕαγ p2ϕαδ p3+w7∑ αβγδµ ϕαγ p1ϕβγ p2ϕδµ p3+w8∑ αβγδµν ϕαβ p1ϕγδ p2ϕµν p3] where ∑′ p1p2p3 means that the sum is restricted to p1+p2+p3= 0. They found a non trivial fixed point below six dimensions. As a test, they recover the previous results of Bray and Roberts. They computed the critical exponent, ν, related with the A and R sectors. However they were unable to
110 CHAPTER 4. MICROSCOPIC DYNAMICS OF THE 3D SPIN GLASS IN PRESENCE OF A MAGNETIC FIELD relate these two νexponents to the physical critical exponents, like γof the spin-glass susceptibility and the standard νof the correlation length. They stated that this identification is difficult since the system presents two mass scales. Some years later [117], Temesv´ari, computed the value of the eight different cubic couplings as a function of the original ones which appear in the Edwards-Anderson Hamiltonian, completing the work started in Ref. [116]. A more recent paper by Bray and Moore [118] states that the de AlmeidaThouless line should disappeared just at six dimensions: h2 AT ∝(6 −D) as D→6.(4.16) They further argue that the break point, x1, of P(q) in the Parisi’s solution should be zero below D < 6. Their final conclusion is that no AlmeidaThouless line can be found below or at six dimensions. Finally, in a recent paper of the Janus collaboration [119], simulations on four dimensional Ising spin glass in presence of an external magnetic field have been developed, showing the presence of a phase transition (please, see Section 7.5 for more details).
CHAPTER 4. MICROSCOPIC DYNAMICS OF THE 3D SPIN GLASS IN PRESENCE OF A MAGNETIC FIELD 111 Dynamics of the D= 3 spin glass in an external magnetic field M. Baity-Jesi, R. Alvarez Ba˜nos, A. Cruz, L. A. Fernandez, J. M. Gil-Narvion, A. Gordillo-Guerrero, D. I˜niguez, A. Maiorano, F. Mantovani, E. Marinari, V. Martin-Mayor, J. Monforte-Garcia, A. Mu˜noz Sudupe, D. Navarro, G. Parisi, S. Perez-Gaviro, M. Pivanti, J. J. Ruiz-Lorenzo, S. F. Schifano, B. Seoane, A. Tarancon, R. Tripiccione and D. Yllanes. To be published 4.2 Introduction The glass transition is a ubiquitous but still mysterious phenomenon in condensed matter physics [120, 121, 122]. Indeed, many materials such as spin glasses, fragile molecular glasses, polymers or colloids display a dramatic increase of relaxation times when cooled down to their glass temperature, Tg. However, the dynamic slowing down is not accompanied by dramatic changes on structural or thermodynamic properties. In spite of this, quite general arguments suggest that the sluggish dynamic must be correlated with an increasing length scale [123]. However, this putative length scale can be fairly difficult to identify. In fact, it was suggested long ago that the slowdown is caused by the collective movements of an increasing number of elements in the system, with a free energy barrier growing with the size of the cooperative regions [124]. These cooperative regions become larger as the temperature gets closer to Tg. The rather recent experimental evidence for cooperative dynamics comes from dynamical heterogeneities [125] or non-linear susceptibilities [126]. The work of Ref. [127] suggests that characteristic length-scales will soon be investigated as well in non-equilibrium, aging materials. It is clear that simple model systems can be a blessing for the study of such a difficult problem. To some extent, spin glasses (which are disordered magnetic alloys [128]) can be such a model system. Upon cooling, they undergo a dynamic slowdown without developing any recognizable magnetic ordering pattern. Their study offers experimental advantages. Time-dependent magnetic fields are a very flexible tool to probe their dynamic response, which can be very accurately measured with a SQUID (for instance, see Ref. [198]). On the theoretical side, they are simple to model, which greatly eases numerical simulation. In fact, special-purpose computers have been built for the simulation of spin glasses [130, 131, 202, 133]. It is then not surprising that the study of spin glasses is ahead in some respects: •We know that the dynamic slowdown is due to a thermodynamic phase
112 CHAPTER 4. MICROSCOPIC DYNAMICS OF THE 3D SPIN GLASS IN PRESENCE OF A MAGNETIC FIELD transition at Tc=Tg[134, 135, 136]. The issue is subtler for supercooled liquids, as we discuss below. •Experiments can measure the size of the glassy magnetic domains, ξ(tw) [137, 138]. These domains are rather large, of the order of 100 lattice spacings [137], compared with any length scale identified for structural glasses [126]. •The Janus dedicated computer [139] allows us to simulate non-equilibrium dynamics from picoseconds to a tenth of a second [133, 199], and to compute equilibrium correlation functions for large lattices and low temperatures [200]. As a result we are able to relate non-equilibrium correlation functions (at finite times) with their equilibrium counterpart in systems of finite sizes [201] (see also Ref. [143]). However, not all is well. We know that spin glasses differ from structural glasses in, at least, two significant respects. First, like all magnetic systems, spin glasses enjoy time-reversal symmetry in the absence of an applied magnetic field. And second, free-energy barriers grow logarithmically with ξ(tw) in spin glasses [199], rather than with a power law as in fragile glasses. The correspondence between spin glasses and structural glasses is more accurate, specially in the mean-field approximation, if one considers instead a rather artificial spin-glass model, the p-spin glass model, with p-body interactions [144, 145]. For odd p, the time-reversal symmetry is broken. The odd-pmodels, at least in the mean-field approximation, display a dynamic phase transition in their paramagnetic phase. Reaching thermal equilibrium becomes impossible in the temperature range Tc< T < Tg. The dynamic transition at Tgis identical to the ideal Mode Coupling transition of supercooled liquids [146]. The thermodynamic phase transition at Tcis analogous to the ideal Kauzmann’s thermodynamic glass-transition [122]. The thermodynamic transition is very peculiar: although it is of the second order (in the Eherenfest sense), the spin-glass order parameter jumps discontinuously at Tcfrom zero to a non-vanishing value. However, the analogy between structural glasses and p-spin glasses was established only in the mean-field approximation. Mean-field is to be trusted only for spatial dimensions larger than the so-called upper critical dimension du. There is no doubt that du>3, hence it is legitimate to wonder how much of the analogy carries out to our three-dimensional world. On the one hand, for supercooled liquids, the ideal Mode Coupling transition is actually a crossover. The power-law divergences predicted by Mode Coupling theory hold when the equilibration time lies in the range 10−13 s< τ < 10−5s. Fitting to those power-laws, one obtains a Mode Coupling temperature TMC.
CHAPTER 4. MICROSCOPIC DYNAMICS OF THE 3D SPIN GLASS IN PRESENCE OF A MAGNETIC FIELD 113 However, τis finite at TMC (typically TMC is a 10% larger than the glass temperature Tgwhere τ∼104seconds). A theory for a thermodynamic glass transition at Tc< Tghas been put forward [147, 148, 149, 150], but it has still not been validated (however, see Ref. [151]). On the other hand, little is known on the behaviour of the p-spin glass model for dimensions below du. A different route to a simple enough model system is quite obvious: break time-reversal symmetry by placing a standard (as opposed to p-spin) spin glass in an external magnetic field. According to mean field [152, 153], though, breaking time reversal is not enough. The mean-field prediction is that, for standard spin glasses on a field, Tc=Tg. Furthermore, the spinglass order parameter would behave continuously when Tcrosses Tc. However, these objections have been challenged for three-dimensional systems (recall that du= 6 [154]). An effective field-theory computation predicts that the spin glass in a magnetic field is the physical realization of a p-spin glass model for spatial dimensions below du[155]. Furthermore, an effective spin glass Hamiltonian in a field has been recently derived for a binary liquid mixture [156]. In fact, whether spin glasses in a magnetic field undergo a phase transition has been a long-debated and still open question (see, e.g., Refs. [157, 158]). Yet, recent numerical simulations in three dimensions [159, 160] did not find the thermodynamic transition predicted by Mean-Field. Experimental studies have been conducted as well, with conflicting conclusions [161, 162, 163, 164]. Only in four dimensions (note that 4 < du= 6) clear signatures of the transition have been found up to now. This exploit required the introduction of special finite-size analysis techniques as well as the power of the Janus special-purpose computer [165]. Our scope here is to explore the dynamical behaviour of three-dimensional spin glasses in a field using the Janus computer. We shall study lattices of size L= 80, where we expect finite-size effects to be negligible [133]. Our time scales will range from 1 picosecond (i.e., one Monte Carlo full lattice sweep [128]) to 0.1 seconds. Hence, if the analogy with structural glasses put forward by Moore and Drossel [155] holds, we should be able of identifying the Mode Coupling crossover. A bonus of studying spin glasses rather than structural glasses is a rather deep theoretical knowledge of the relevant correlation functions [166]. Hence, we shall be able to correlate the equilibration time τwith the correlation length ξ.
114 CHAPTER 4. MICROSCOPIC DYNAMICS OF THE 3D SPIN GLASS IN PRESENCE OF A MAGNETIC FIELD 4.3 Model and observables 4.3.1 Model We studied a three dimensional cubic lattice system with volume V=L3(L being the linear size) and periodic boundary conditions. On every node of the lattice there is an Ising spin, σx=±1 and nearest neighbors are joined bye quenched bimodal couplings, Jxy =±1. We also include a local magnetic field, hx, on every node. The magnetic field is Gaussian distributed with zero mean and variance H. Instead of continuous values, we used discrete values for the field by using the Hermite integrals of its probability distribution [167] (see the Appendix for more details on the implementation). We made this transformation to use more efficiently the supercomputer Janus [139, 202, 203]. We checked the compatibility of our approach by comparing with real Gaussian fields simulated on PCs (see also the Appendix). The Hamiltonian of the model is H=−∑ ⟨x,y⟩ Jxy σxσy−∑ x hxσx,(4.17) where ⟨x,y⟩means sum over nearest neighbors. A given realization of couplings, Jxy , and external field, hx, defines a sample. We have simulated four replicas in parallel with the same couplings and protocols (annealing and direct quench). 4.3.2 Observables First, a couple of useful definitions of local quantities. On every node xof the lattice we have the local overlap: qx(t) = σ(1) x(t)σ(2) x(t),(4.18) where the superscripts are the replica indices. The total overlap is written as q(tw) = 1 V∑ x qx(tw),(4.19) where (···) means sample average (over the J’s and h’s). Notice that lim tw→∞q(tw) = qmin ,(4.20) where qmin is the minimum overlap allowed by the system. In addition, we have focused in this work on the magnetic energy defined as Emag(tw) = 1 V∑ x hxσx(tw).(4.21)
CHAPTER 4. MICROSCOPIC DYNAMICS OF THE 3D SPIN GLASS IN PRESENCE OF A MAGNETIC FIELD 121 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1001021041061081010 1012 q(t) and W=1-T*Emag/H2 t q(t) W Figure 4.9: q(tw) and W(tw) for H= 0.3 and base = 105(annealing run). Every step in the figure corresponds with a change of the temperature. 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0 0.005 0.01 0.015 0.02 0.025 0.03 0.035 0.04 W(t)-q(t) t-0.22 H=0.1 H=0.2 H=0.3 Figure 4.10: Extrapolation of the difference between W(t) and q(t) as a function of a power of time, for the three external magnetic fields simulated. Bottom to top: H= 0.1, 0.2 and 0.3. Temperature T= 0.5. for the three magnetic fields, hence this temperature at these three magnetic fields behaves as driven by a spin glass phase. The dependence of the data plotted in these figures with the base parameter of the annealing procedure
122 CHAPTER 4. MICROSCOPIC DYNAMICS OF THE 3D SPIN GLASS IN PRESENCE OF A MAGNETIC FIELD 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0 0.005 0.01 0.015 0.02 0.025 0.03 0.035 0.04 0.045 W(t)-q(t) t-0.22 H=0.1 H=0.2 H=0.3 Figure 4.11: Extrapolation of the difference between W(t) and q(t) as a function of a power of time, for the three external magnetic fields simulated. Bottom to top: H= 0.1, 0.2 and 0.3. Temperature T= 0.6. is inside the error bars. If we examine the next higher temperature, T= 0.6, we can observe that the data corresponding to H= 0.3 extrapolate to zero, while the two lower magnetic fields, H= 0.1 and 0.2, still have a non-zero extrapolated value: we can conclude that the point (T, H) = (0.6,0.3) is just in the paramagnetic phase, whereas the pairs (0.6,0.2) and (0.6,0.1) are still in a spin glass phase. In particular, we can state that the spin glass phase (the de AlmeidaThouless line) satisfies T1(H= 0.3) >0.5. For T= 0.8 (see Fig. (4.13)) only H= 0.1 extrapolates to a non-zero value, whereas at T= 0.9 (see Fig (4.14)) all three magnetic field extrapolates to a non positive value. One can roughly estimate that T1(H= 0.3) ≃0.6, T1(H= 0.2) ≃0.7 and T1(H= 0.1) ≃0.8. We can study in more detail the dependence of the power law exponent with the temperature. As it has been described above we have fitted the difference between q(t) and W(t) following the power law described by Eq. (4.28). This equation, with a≥0 should hold only deeply in the spin glass phase. If we approach the transition from below we would start to see the critical effects of the (thermodynamical) critical point, and the exponent x begin to be controlled by this critical point and not by the “critical” spin glass phase (Goldstone phase). So, in the critical region we should expect: W(t)−q(t) = f txc,(4.31)
CHAPTER 4. MICROSCOPIC DYNAMICS OF THE 3D SPIN GLASS IN PRESENCE OF A MAGNETIC FIELD 123 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0 0.01 0.02 0.03 0.04 0.05 0.06 W(t)-q(t) t-0.22 H=0.1 H=0.2 H=0.3 Figure 4.12: Extrapolation of the difference between W(t) and q(t) as a function of a power of time, for the three external magnetic fields simulated. Bottom to top: H= 0.1, 0.2 and 0.3. Temperature T= 0.7. 0 0.05 0.1 0.15 0.2 0.25 0.3 0 0.01 0.02 0.03 0.04 0.05 0.06 W(t)-q(t) t-0.22 H=0.1 H=0.2 H=0.3 Figure 4.13: Extrapolation of the difference between W(t) and q(t) as a function of a power of time, for the three external magnetic fields simulated. Bottom to top: H= 0.1, 0.2 and 0.3. Temperature T= 0.8. where in general xc(driven by the critical point) should be different from x (driven by the spin glass phase, which is a critical one). Finally well above
124 CHAPTER 4. MICROSCOPIC DYNAMICS OF THE 3D SPIN GLASS IN PRESENCE OF A MAGNETIC FIELD 0 0.05 0.1 0.15 0.2 0.25 0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 W(t)-q(t) t-0.22 H=0.1 H=0.2 H=0.3 Figure 4.14: Extrapolation of the difference between W(t) and q(t) as a function of a power of time, for the three external magnetic fields simulated. Bottom to top: H= 0.1, 0.2 and 0.3. Temperature T= 0.9. the critical region we should expect a stretched exponential behavior (see Eq. (4.29)). From the previous discussion, and assuming the onset of a phase transition, it is clear that the x-exponent should take a constant value at lower temperatures (here we are assuming that the phase transition is Universal in the magnetic field), then change as we reach the critical region, and finally change again in the high temperature region since the pure power law is not longer valid (the behavior should switch to a stretched exponential). If we try to fit the high-temperature region with Eq. (4.28) with a= 0 or Eq. (4.31) instead of, for example Eq. (4.29), we will obtain a higher value of the x-exponent to compensate the lack of the exponential. This is just what happens in Fig. (4.15). As an additional test we can monitor the behavior of the constant term in the power law fit (see Eq. (4.28). We present the dependence of awith temperature in Fig. (4.16), and we can observe that above a threshold (which depends on the magnetic field) the value of astarts to be negative. This conclusion reinforces the results obtained using a constant x-value in the fits. Notice that the values of xand apresented in figures (4.15) and (4.16) are obtained doing a three-parameter fit on our data, and that our choice for x (which is 0.22) is compatible with the xfrom the three-parameter fit in the low temperature region and for the three values of the magnetic fields.
CHAPTER 4. MICROSCOPIC DYNAMICS OF THE 3D SPIN GLASS IN PRESENCE OF A MAGNETIC FIELD 125 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 0.55 0.6 0.65 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 x T H=0.30 H=0.20 H=0.10 Figure 4.15: Exponent of the extrapolation of the difference between W(t) and q(t), x(see Eq. (4.28)), as a function of temperature, for the external magnetic fields simulated. Value computed from a three-parameter fit. -0.02 0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 a T H=0.30 H=0.20 H=0.10 Figure 4.16: Asymptotic value of the extrapolation of the difference between W(t) and q(t), a(see Eq. (4.28)), as a function of temperature, for the external magnetic fields simulated. Value computed from a three-parameter fit. We can use the temperature at which abecomes negative as our estimate of T1(H). From figure 4.16 we can estimate, now leaving vary the exponent
126 CHAPTER 4. MICROSCOPIC DYNAMICS OF THE 3D SPIN GLASS IN PRESENCE OF A MAGNETIC FIELD x,T1(H= 0.3) ≃0.65(5), T1(H= 0.2) ≃0.80(5), T1(H= 0.1) ≃0.96(5). It is clear that we have found a spin glass region in field, nevertheless the method used cannot allow us to obtain a precise value of the de AlmeidaThouless line. 4.5.2 High temperature Region: Computation of the relaxation times Once we have some estimates of T1(H) we can try to study the dynamical behavior of the system in the high temperature region. We have observed that our data for W(tw)−q(tw) in this high temperature region follows very well the stretched exponential behavior (see Eq. (4.29)). This allows us to compute, using our annealing runs, the relaxation time as a function of the temperature. We should keep in mind that the computed time should be less than the maximum time the system is in a given temperature during the annealing process. 1 100000 1e+10 1e+15 1e+20 1e+25 0.6 0.8 1 1.2 1.4 1.6 1.8 2 τ T H=0.3 H=0.2 H=0.1 Annealing max. times Figure 4.17: Behavior of the correlation time (τ) as a function of the temperature for the three magnetic fields simulated (see Eq. (4.30)). We also plot the best fits we have had using the critical law of τ, see the text for more details. Finally we have a discountinous line of triangles which marks the maximum times simulated during the annealing procedure at a given temperature, which marks a cutoff on our computation of the relaxation times. Computing a fit to Eq. (4.29) is difficult due to the extreme correlation of our data, which prevents us from inverting its full covariance matrix (neces-
CHAPTER 4. MICROSCOPIC DYNAMICS OF THE 3D SPIN GLASS IN PRESENCE OF A MAGNETIC FIELD 127 sary to define the χ2goodness-of-fit indicator). Therefore, we consider only the diagonal part of the matrix in order to minimize χ2and take correlations into account by repeating this procedure for each jackknife block in order to estimate the errors in the parameters. This is, of course, only an empirical procedure, but one that has been shown to work well under these circumstances (see, e.g., Ref. [199], especially sections 2.4 and 3.2). In Fig. (4.17) we show the computed relaxation time as a function of the temperature and for the three simulated magnetic fields. In addition, we have plotted the maximum times the system spends at each temperatures (for our largest value of base). 0.3 0.5 0.7 0.9 0.8 1 1.2 1.4 1.6 1.8 2 β T H = 0.1 H = 0.2 H = 0.3 Figure 4.18: Behavior of the stretching exponent β(see Eq.(4.29)) as a function of Tfor our three simulated magnetic fields. Let us mention, finally, that a possible additional source of uncertainty in our determination of τis the depedence of the fit on the value of β. Indeed, for each Twe are fitting simultaneously for x,A,τand βin (4.29). However, a small variation in βcan have a large effect on τ, which may lead us to think that the fit is unstable and unreliable. Fortunately (see Fig. 4.18), β is actually a very smooth monotonic function of T, which leads us to believe that our determination of τis sound. Fig. (4.17) shows us that the relaxation times are diverging very quickly, as a function of temperature, for the three magnetic field. One can examine in detail if this behavior is driven by a divergence at finite temperature. Again, and following Ref. ([131]) we try to fit the “divergence” of the relaxation time under the onset of a phase transition at finite temperature: e.g. by
128 CHAPTER 4. MICROSCOPIC DYNAMICS OF THE 3D SPIN GLASS IN PRESENCE OF A MAGNETIC FIELD using Eq. (4.30). The continuous lines in the fit correspond to this kind of fit (with very good χ2/dof, dof being the number of degrees of freedom). We have obtained the following values: •H= 0.1: T2= 0.97(6) and zν = 5.8(7). Using only T≥1.2 [χ2/dof = 0.77]. •H= 0.2: T2= 0.72(6) and zν = 7.3(1.0). Using only 0.9≤T≤1.7 [χ2/dof = 0.79]. •H= 0.3: T2= 0.66(8) and zν = 6.2(1.6). Using only 0.8≤T≤1.7 [χ2/dof = 0.5]. Since we have computed the correlation length in addition to the relaxation times, we can address the issue of the dependence of τwith ξ, which gives us useful information on the dynamics. 10 100 1000 10000 100000 1e+06 1e+07 1e+08 1e+09 1e+10 1e+11 1 10 τ ξ12 H=0.3 H=0.2 H=0.1 Slope H=0 Figure 4.19: Behavior of τagainst the correlation length (ξ12) for the three magnetic fields simulated. We have also marked the H= 0 behavior: τ≃ξz, with z= 6.86. In Fig. (4.19) we plot τagainst ξ12 (defined using Eq. (4.26)) for the three simulated magnetic fields. In addition, to control, we have plotted the H= 0 behavior (τ∝ξz, using the critical temperature value for the dynamical critical exponent z≃6.9). From this figure, it is clear that the dependece of τwith ξhave changed when we have turned on the external magnetic field, however, we lack of accuracy in order to determine the analytical dependence.
CHAPTER 4. MICROSCOPIC DYNAMICS OF THE 3D SPIN GLASS IN PRESENCE OF A MAGNETIC FIELD 129 4.6 Discussion of the Results Since in our high-temperature study we have followed closely Ogielski’s approach, we need to put the exponents and critical temperatures computed by Ogileski (remember at H= 0) in relation with the most accurate values found in the literature, in order to asses our own data. •Ogileski provided as a critical temperature Tg= 1.175(25) which should be compare with that computed in Ref. [170]: Tc= 1.109(10): monitoring the relaxation times gives us an overestimated value (6%) of the critical temperature. •He obtained zν = 7.0(8) and 7.9(1). The most recent and accurate values for ν= 2.53(8) [170] and z= 6.86(16)[199], providing us zν = 17.4(7). Hence, the computed value of zν is off by a factor of two. •Experimentalists have also followed this strategy since they are able to measure q(t) in the high-temperature region. The experimental value for Ising spin glass can be quoted as: zν ≃10.5(1.0) [169]. Summarizing, (at H= 0) one obtains an overestimated temperature (+6%) and a factor two off value for the product νz. Notice the robustness of the procedure even providing a wrong value of the product, essentially the same number (near 10) is obtained in real experiments. It is clear that the lack of corrections-to-scaling in the analysis of the relaxation times has strong effect in zν but not so much in the critical temperature. With the available computational facilities we are unable to improve this procedure. Once we have discussed the methodology (and some drawbacks) used to obtain our data, we can try to put a coherent physical picture. At this point, we have clearly three possibles scenarios: NPT. No phase transition at all. Both T1(H) and T2(H) should eventually drift to zero temperature or be the effect from a crossover from the H= 0 phase transition. The experiments in a field should give T2(H) (they can not access T1(H) in this way). They obtain a good dynamical scaling in field with zν ≃9, but this is usually interpreted as a crossover effect. In addition, the low-temperature behavior of the difference (1/t0.22) could change if we simulate long waiting times (however, we remind the reader that we are already simulating up to the beginning of the experimental time scales). In addition, the value of T1(H= 0.1) is compatible with the critical temperature of the model n absence of magnetic field (Tc(H= 0) ≃1.1,
130 CHAPTER 4. MICROSCOPIC DYNAMICS OF THE 3D SPIN GLASS IN PRESENCE OF A MAGNETIC FIELD so our data for H= 0.1 coild be strongly affected by the H= 0 critical point. Yet, the values for T1(H= 0.2) and T1(H= 0.2) are not near to T≃1.1, so, in principle, these two magnetic field should be avoided the crossover effect of the H= 0 critical point. 2PT. T1(H)< T2(H) Scenario. This scenario is the most suggestive one from the point of view of the hypothetical correspondence beteween structural glasses and spin glasses on a magnetic field [155, 156]. The replica theory for structural glasses [147, 148, 149, 150] suggests that T2(H) would rather correspond to the Mode Coupling temperature (which is rather a crossover in three spatial dimensions), while a real thermodynamic phase transition would take place at T1(H). We note that the existence of a thermodynamic glass transition is being vigorously debated by the supercooled liquids community [122]. 1PT. Only one thermodynamical phase transition. In this light, T1(H) = T2(H), since the phase transition drives the divergence of the relaxation times and also the breakdown of the law given by Eq. (4.28). The main difference of this work regarding dynamical experimental studies is that we can also compute T1(H) in addition to T2(H). We note that this scenario is the one predicted by Mean Field theory. Regarding the last two scenarios, our values of T1(H) and T2(H) are very similar, but they are not so accurate to fix the possible escenario. As cited in the fisrt scenario, we cannot even discard an eventual crossover of both T1 and T2to zero (simulating larger values of base). 4.7 Conclusions We have tried to characterize the behavior of the three-dimensional spin glass both in the high and low temperature regions, monitoring the behavior of the difference W(t)−q(t). These studies have allowed us to determine (at least in our range of simulated times) two changes of regime as a function of temperature: in the first one (T1(H)), below that, the low temperature phase behaves as a spin glass one; in the second one (T2(H)), we have obtained a divergence of the relaxation times. Numerically we have found that Ts(H) is roughly similar to T1(H). We are simulating the beginning of the experimental times, so, we will have the same advantages and drawbacks as in real experiments: in particular, we cannot discard a change of the low-temperature exponent (x≃0.22), which, eventually, can drive T1(H) to zero.
List of Figures 1.1 Susceptibilidad del CuMn con un protocolo de enfriamiento en campo magn´etico y otro de enfriamiento en campo cero. Figura de la Ref. [10]. . . . . . . . . . . . . . . . . . . . . . . 3 1.2 Magnetizaci´on remanente del (Fe0.15Ni0.85)75P16B6Al3. Figura delaRef.[11]. .......................... 3 1.3 Representaci´on esquem´atica de la distribuci´on del overlap en la fase paramagn´etica. . . . . . . . . . . . . . . . . . . . . . . 9 1.4 Representaci´on esquem´atica de la soluci´on hallada para RSB. . 13 1.5 Representaci´on esquem´atica de la distribuci´on del overlap de lasoluci´onRSB........................... 13 1.6 Representaci´on esquem´atica de la distribuci´on del overlap en un droplet.............................. 15 1.7 Un ejemplo de una plaqueta 2 ×2 frustrada. . . . . . . . . . . 19 1.1 Susceptibility of CuMn with a field-cooling protocol and a zero-field cooling. Figure from Ref. [10]. . . . . . . . . . . . . 23 1.2 Remanent magnetization of (Fe0.15Ni0.85)75P16B6Al3. Figure fromRef.[11]............................ 23 1.3 Schematic representation of the distribution of the overlap in the paramagnetic phase. . . . . . . . . . . . . . . . . . . . . . 28 1.4 Schematic representation of the solution found to RSB. . . . . 33 1.5 Schematic representation of the distribution of the overlap the RSBsolution............................ 33 1.6 Schematic representation of the distribution of the overlap in adroplet. ............................. 35 1.7 An example of a 2 ×2 frustrated plaquette. . . . . . . . . . . 38 2.1 Schematic representation of the solution found for p < p∗. . . . 43 2.2 Schematic representation of the solution found for p > p∗and T2< T < Tc............................. 45 233
234 LIST OF FIGURES 2.3 Schematic representation of the solution found for p > p∗and T < T2< Tc............................. 46 2.4 Evolution of the viscosity of several liquids. Notice that, although the evolution is different, all of them reach the same value of the viscosity. This figure is the famous Angell plot, fromRef.[45]........................... 47 2.5 Log-binning thermalization test for p= 5. For all data points the point size is bigger than the corresponding error bar. . . . 57 2.6 As in figure 2.5, but p=6..................... 58 2.7 The autocorrelation function (2.31) for one generic sample (p= 6, L=8). .......................... 59 2.8 Integrated autocorrelation time, τint, for all p= 5, L= 8 samples. τint is in units of blocks of ten measurements, i.e. of 20103 MCS. Samples above the green line have been “extended” (see the text for a discussion of this issue). . . . . . . . . . . . . . . 60 2.9 Overlap correlation length in lattice size units as a function of the inverse temperature βfor L= 4, 6, 8 and 12. Here p= 5. . 61 2.10 As in figure 2.9, but p=6..................... 62 2.11 Magnetic susceptibility as a function of βfor L= 4, 6, 8 and 12. Here p=5. .......................... 64 2.12 As in figure 2.11, but p=6. ................... 65 2.13 In the bottom plot: βcversus p, and the straight line f(p) = p. Middle plot: νas a function of p. We also show (dashed line) the value which marks the onset of a disordered first order phase transition (νfirst = 2/3). Upper plot: ηqas a function of p. 67 3.1 (Color online) The quantity XTas defined in the text, as a function of qfor lattice size L= 24 (top) and L= 32 (bottom) at temperature T≃0.64Tc. Insets show a magnified view of the region q∼0.6 (log-log plot). Plots show data for XTcomputed only with triplets of independent configurations (ABC), with triplets in which two configurations belong to the same Monte Carlo history (AAB), and triplets in which all configurations come from the same Monte Carlo history (AAA). No significant difference shows up as long as we take enough uncorrelated configurations from the same replica. . . 94
LIST OF FIGURES 235 3.2 (Color online) Top: X2as a function of the corresponding polynomial in X1(Eq. (3.76)). The straight line is the theoretical prediction (unit slope). Center: the ratio X2/X1as a function of X1, where the straight line is the theoretical prediction. Bottom: the squared difference K2= [X2−(X1+ 2X2 1)/3]2as function of X1. Data refer to T∼0.64Tc............ 95 3.3 (Color online) Data at T∼0.64Tc. Top: X3as a function of the corresponding polynomial in X1and XT(Eq. (3.77)). The straight line is the theoretical prediction (unit slope). Bottom: the squared difference K3= [X3−(2XT+ 2X1+ 6X2 1+ 5X3 1)/15]2 as function of X1,T= 0.64Tc. Lines connecting points are only a guide to the eye. . . . . . . . . . . . . . . . . . . . . . . 96 3.4 (Color online) Top: The squared difference [XT−X2 1]2as a function of X1. Bottom: the quantity Ku 3= [X3−(2X1+ 8X2 1+ 5X3 1)/15]2 as a function of X1. All data for T∼0.64Tcand for lattice sizes L= 16,24,32. The lines connecting the data points are only intended as a guide to the eye. . . . . . . . . . . . . . . . 97 3.5 (Color online) Square difference [XT−X2 1]2(left) and the quantity Ku 3= [X3−(2X1+ 8X2 1+ 5X3 1)/15]2(right) as a function of X1. Top: for T= 0.75Tcand L= 8,16,24,32. Bottom: for T= 0.57Tcand L= 8,16,24. .................. 98 3.6 (Color online) Top: The conditioned probability P(q12|q34) (open squares) for L= 32 and T∼0.64Tcand two values of q34 = 0.211 (left) and q34 = 0.367 (right). We also plot 2P(q12)/3 (open circles) and the difference (full triangles) of the two above quantities (Eq. (3.85) in the text), scaled by a factor 2 for a better view. q34 and qEA values are indicated by vertical lines for visual reference. We took the value qEA(L= 32, T = 0.64Tc)∼0.72 as given in Ref. [89]. Bottom: The difference P(q12|q34)−2P(q12)/3 with q34 = 0.367, for different lattice size compared at temperatures T= 0.75Tc,T= 0.64Tc, T= 0.57Tc. ............................ 99 3.7 (Color online) Asymptotic behavior of the cumulative probability ΠC q(z) (Eq. (3.81)). Top: small-zdecay for L= 32, T= 0.64Tcand q= 0.3125. Bottom: comparison of the exponent x(q) obtained by the two methods described in the text (uppermost data points represent values obtained by fitting ΠC q(z), lowermost data points come from integrating the P(q)), for some lattice sizes, many cut-off values qand temperatures T∼0.57Tc(left) and T∼0.64Tc(right). . . . . . . 100
236 LIST OF FIGURES 3.8 (Color online) Comparison between the Monte Carlo data of the P(q) and the convolution computed as described in the text (solid lines). Top: L= 32, T∼0.64Tcand T∼0.75Tc. Center: L= 24, T∼0.57Tcand T∼0.64Tc. Bottom: the conditioned probability P(q12|q34 =q0) for L= 32, T∼0.64Tc and some values of q0. ......................101 4.1 Phase diagram in T−hvariables in the RSB scenario. The de Almeida-Thouless line separates the paramagnetic and spin glassphases.............................104 4.2 Probability distribution of qin presence of an external magnetic field in RSB scenario. . . . . . . . . . . . . . . . . . . . . 104 4.3 Probability distribution of qin presence of an external magnetic field in Droplet scenario. . . . . . . . . . . . . . . . . . . 105 4.4 Behavior for the dynamical overlap, q(τ) (which is proportional to the quantity plotted in the y-axis), as a function of time for different magnetic fields. Figure from Ref. [108]. . . . 106 4.5 Relative decrease of Tc(h)/Tc(0) with increase field for ρ= 1.5 and h= 0, 0.1, 0.15 and 0.2 versus the relative decrease of χ∗ (ZFC susceptibility). Figure from Ref. [110]. Experimental data from Fe0.5Mn0.5TiO3, see Ref. [108]. . . . . . . . . . . . . 106 4.6 q(tw) and W(tw) at T= 0.7 and H= 0.1.............119 4.7 q(tw) and W(tw) at T= 0.7 and H= 0.3. ...........119 4.8 q(tw) and W(tw) for H= 0.1 and base = 105(annealing run). Notice that every step in the figure corresponds with a change ofthetemperature. .......................120 4.9 q(tw) and W(tw) for H= 0.3 and base = 105(annealing run). Every step in the figure corresponds with a change of the temperature...............................121 4.10 Extrapolation of the difference between W(t) and q(t) as a function of a power of time, for the three external magnetic fields simulated. Bottom to top: H= 0.1, 0.2 and 0.3. Temperature T= 0.5..........................121 4.11 Extrapolation of the difference between W(t) and q(t) as a function of a power of time, for the three external magnetic fields simulated. Bottom to top: H= 0.1, 0.2 and 0.3. Temperature T= 0.6..........................122 4.12 Extrapolation of the difference between W(t) and q(t) as a function of a power of time, for the three external magnetic fields simulated. Bottom to top: H= 0.1, 0.2 and 0.3. Temperature T= 0.7..........................123
LIST OF FIGURES 237 4.13 Extrapolation of the difference between W(t) and q(t) as a function of a power of time, for the three external magnetic fields simulated. Bottom to top: H= 0.1, 0.2 and 0.3. Temperature T= 0.8..........................123 4.14 Extrapolation of the difference between W(t) and q(t) as a function of a power of time, for the three external magnetic fields simulated. Bottom to top: H= 0.1, 0.2 and 0.3. Temperature T= 0.9..........................124 4.15 Exponent of the extrapolation of the difference between W(t) and q(t), x(see Eq. (4.28)), as a function of temperature, for the external magnetic fields simulated. Value computed from a three-parameter fit. . . . . . . . . . . . . . . . . . . . . . . . 125 4.16 Asymptotic value of the extrapolation of the difference between W(t) and q(t), a(see Eq. (4.28)), as a function of temperature, for the external magnetic fields simulated. Value computed from a three-parameter fit. . . . . . . . . . . . . . . 125 4.17 Behavior of the correlation time (τ) as a function of the temperature for the three magnetic fields simulated (see Eq. (4.30)). We also plot the best fits we have had using the critical law of τ, see the text for more details. Finally we have a discountinous line of triangles which marks the maximum times simulated during the annealing procedure at a given temperature, which marks a cutoff on our computation of the relaxation times.126 4.18 Behavior of the stretching exponent β(see Eq.(4.29)) as a function of Tfor our three simulated magnetic fields. . . . . . 127 4.19 Behavior of τagainst the correlation length (ξ12) for the three magnetic fields simulated. We have also marked the H= 0 behavior: τ≃ξz, with z= 6.86..................128 4.20 Thermal energy (E), magnetic energy (W(t) and overlap (q(t)) as a function of time for L= 8, T= 0.7 and H= 0.3. We have plotted the results from a fully Gaussian (G.), n= 2 and n= 5 numerical simulations. Notice that all three simulations provided us with the same values of these three observables. . 132 5.1 Locus of the Fisher zeros of a 2D Ising model, L= 16. Figure fromRef.[177] ..........................134
238 LIST OF FIGURES 5.2 The four first zeros at β=βc. In order to appreciate the scaling better, we show only the data for L≥12 and compare to equation (5.12), fixing x2=ω= 1.0(1) from [212] and performing a global fit for a common value of x1(see text). We obtain x1= 2.67(6)[1], with a chi-square per degree of freedom of χ2/d.o.f.= 5.88/7. .................145 5.3 Scaling of the zeros at β= 1.4, with a best fit to (5.25) for L≥16. We obtain x1= 2.842(11), with χ2/d.o.f.= 7.34/7. . . 150 5.4 Scaling of the zeros at β= 1.2, with a best fit to (5.25) for L≥16. We obtain x1= 2.844(10), with χ2/d.o.f.= 2.89/7. . . 150 5.5 χ/V =⟨q2⟩versus the lattice size for β= 1.2 and 1.4. Notice that none of the temperatures have reached the plateau asymptoticvalue..........................151 5.6 Integrated density of zeros versus the zeros at the critical point. a2= 1.16(2). .......................151 5.7 Integrated density of zeros versus the zeros for β= 1.2. . . . . 152 5.8 Integrated density of zeros versus the zeros at β= 1.4. . . . . 152 5.9 Integrated density of the zeros, for the largest lattice L= 32 and the lowest temperature β= 1.4. Notice that we are almost, but not in, the linear regime. The data are well fitted with b= 1.068(10). ........................153 5.10 Histogram (N(ϵ) versus ϵ) for the 1000 first zeros computed for L= 32 and β= 1.4. Notice the lack of symmetry of the histogram and the presence of events for large values of the zeros.................................154 5.11 Integrated density of the zeros, for the largest lattice L= 32 and temperature β= 1.2 using the average of zeros. We have also plotted the median values. We have marked the expected slope at the origin, using the Edwards-Anderson order parameter computed in Ref. [200] for the L= 32 lattice. . . . . . . . 155 5.12 Integrated density of the zeros, for the largest lattice L= 32 and lowest temperature β= 1.4 using the average of zeros. We have also plotted the median values. We have marked the expected slope at the origin, using the Edwards-Anderson order parameter computed in Ref. [200] for the L= 32 lattice. 156
LIST OF FIGURES 239 6.1 Susceptibility versus the temperature. Solid line is the reference one (without any stop). Open diamonds mark measurements while decreasing temperature with a stop at 12 K during 7h. Solid circles mark measurements increasing the temperature. The rate of the change of the temperature is 0.1 K/min. Figure from reference [219]. . . . . . . . . . . . . . . . 158 6.2 Susceptibility versus the temperature. In this figure, the stops are of 30 min. Figure from Ref. [220]. . . . . . . . . . . . . . . 159 6.3 Susceptibility at t0= 624 and maximum twvs temperature. . 160 6.4 Susceptibility at t0= 390624 and maximum twvs temperature. 161 6.5 Coherence length versus t. ....................161 6.6 Coherence length (at the largest twin every temperature) versustemperature. .........................162 6.7 Coherence length versus t. Simulations of the same samples at fixed temperatures T= 0.9 and T= 0.8 are also plotted. . . 162 6.8 Coherence length versus t. Simulations of the same samples at fixed temperatures T= 0.9 and T= 0.8 are also plotted. . . 163 7.1 Probability distribution of the overlap at temperature T=0.625 Figure from J. Stat. Mech. P06026 (2010) [230]. . . . . . . . . 168 7.2 Probability distribution of the overlap at temperature T=0.703 Figure from J. Stat. Mech. P06026 (2010) [230]. . . . . . . . . 168 7.3 Crossovers of Fq/Lyfor a couple of values of y:y= 2.35 (top) and y= 2 (bottom). The insets show in detail the crossing regions. Figure from Phys. Rev. Lett. 105 177202. . . . . . . 170 7.4 Top: plot of the ξ2correlation length versus temperature at h= 0.15. Any intersection is found. Bottom: plot of R12 versus temperature at h= 0.15. One can observe now intersections. Figure from PNAS 109 6452-6456. . . . . . . . . . . 171 A.1 Configuration of a Janus board . . . . . . . . . . . . . . . . . 181 A.2 AJanusboard...........................182 A.3 Nearest-neighbour toroidal network of the SPs of a Janus board182 A.4 Update of a white spin. One only needs black spins . . . . . . 183 E.1 Evolution of the index of the βwhere a given configuration stays. Notice that . Data from Potts (Section 2) simulations: a sample with p= 5 and L=12. ................202
240 LIST OF FIGURES
List of Tables 2.1 Details of the simulations for p=5................ 55 2.2 Details of the simulations for p=6................ 55 2.3 Numerical values of our estimates for the crossing point of the curves ξ/L. We give βcross, the thermal critical exponent ν, the anomalous dimension of the overlap ηq, and the anomalous dimension of the magnetization ηm................ 62 2.4 As in table 2.3, but p=6. .................... 63 2.5 Critical parameters as a function of p. All data are for binary couplings, with zero expectation value. By Rwe denote the ratio between the critical βin three dimensions and that computed in Mean Field. . . . . . . . . . . . . . . . . . . . . 66 3.1 A summary of parameters of the simulations we have used in this work. For each lattice size, L, we considered NSsamples, with four independent real replicas per sample. For the Parallel Tempering algorithm, NTtemperatures were used between Tmin and Tmax, uniformly distributed in that range (except in the case of L= 8, in which we have 7 temperatures uniformly distributed between 0.435 and 1.575 plus the 3 temperatures 0.150, 0.245 and 0.340). Our MCS consisted of 10 Heat-Bath sweeps followed by 1 Parallel Tempering update. More detailed information regarding these simulations can be found inRef.[89]. ............................ 82 3.2 Temperature values for each lattice size (Tc= 1.109 [104, 105]). 87 3.3 Results of the fitting procedure of Eq. (3.94) on numerical P(q) data, with kernel exponent k= 2.5 (see Eq. (3.91)). All errors on parameters are jackknife estimates. We used the symbol χ2 in the table to denote the sum of squares of residuals, which is not a true chi-square estimator as the values of P(q) at different qare mutually correlated. . . . . . . . . . . . . . . . 91 241
242 LIST OF TABLES 4.1 Details of the simulations at fixed temperature. MCS means total Monte Carlo steps and N means the number of samples simulated. ............................117 4.2 Details of the simulations with the annealing algorithm. The same notation as in Table (4.1) and Tinit and Tend mark the initial and final temperatures of the annealing procedure. . . . 118 5.1 Summary of the simulations. NTis the number of simulated temperatures (evenly spaced between Tmin and Tmax); Nmes is the number of Monte Carlo steps (updates of the whole lattice) between measurements; Nmed HB is the average simulation time (since we use the random-walk technique the simulation time depends on the sample); Nsam is the number of simulated samples. We have simulated four real replicas for each sample. Finally, L= 8 and L= 12 have been simulated on PCs and L= 16, L= 24 and L= 32 on Janus. . . . . . . . . . . . . . . 142 5.2 Fits of the zeros to ϵj(L) = AjL−x1, for L≥Lmin. As we can see, with Lmin = 8 the χ2per degree of freedom is acceptable only for j= 1,2, but with Lmin = 12 all the zeros have a reasonable fit. However, the value of x1grows with j, an indication that we have to consider corrections to scaling (see text). ...............................144 5.3 Scaling of the zeros in the low-temperature phase. For the two considered temperatures (β= 1.2,1.4) we first show a fit without corrections to scaling for L≥16, that is ϵj(L)≃ AjL−x1. As explained in Section 5.5.1, this is a global fit for the four zeros, considering their full covariance matrix. We then consider the same fit with corrections to scaling, trying different values for ω(see the text for more details). In all cases x1is smaller than the expected value x1=D= 3. . . . 146 7.1 Simulation details. . . . . . . . . . . . . . . . . . . . . . . . . 167 7.2 Critical parameters for different values of external magnetic fields.................................171