Enrique Gaztañaga

Inside a Black Hole: the illusion of a Big Bang

Keywords: Cosmology: dark energy, cosmic background radiation, cosmological parameters, early Universe, inflation

whether they are published or not. The documents may come

Inside a Black Hole: the illusion of a Big Bang

INTRODUCTION

A Schwarzschild BH metric (BH.SW) represents a singular object of mass . The event horizon

≡ 2 prevent us from interacting with the inside which makes BHs good Dark Matter candidates. Physically, a singular point does not make sense. 1 But objects with mass and sizes matching have been observed. What is the metric inside? What happens when they accrete matter or when two BHs merge? Do BHs grow and co-evolve with galaxies (e.g. Kormendy & Ho 2013)? Do observed BH form in stellar collapse or are they seeded by primordial BHs? How do primordial BH form (e.g. Kusenko & et al. 2020)? Most of these modelings assume the BH.SW solution, but can we actually answer any of this if we do not have a physical model for the BH interior?

Here, we look for an alternative solution to the BH.SW interior, defined as a non singular classical object of size which reproduces the BH.SW metric for the outside >

. A physical BH of size <

and mass , has a density:

= = 3 -2 8 = 3 -2 32 3 . (1)
This is more compact than any form of regular matter (Buchdahl 1959). The highest known density for a stellar object is that of a Neutron star, which has the density of an atomic nucleus, but is still a few times larger than . Yet BH of different have been observed, even of galactic size. To achieve such a high density for ★ E-mail:gaztanaga@gmail.com 1 This is why it took Newton over 20 years to publish the inverse square law of gravity and re-invent integration on the way (darkcosmos.com). a perfect fluid, the radial pressure inside a BH needs to be negative (Brustein & Medved 2019 and references therein). Cosmologist are used to this type of fluids, which are called Quintessence, Inflation or Dark Energy (DE). So, could the inside of a BH be DE? Mazur & Mottola (2015) have argued that the same DE repulsive force that causes cosmic acceleration could also prevent the BH collapse, resulting in the so call gravastar solution. The simplest DE is the ground state ≡ 0 ( ) of a scalar field ( ). Gaztañaga (2021b) have argued that a constant vacuum energy does not gravitate. But, as we will show here, a false vacuum (FV) discontinuity = 0 +Δ does gravitate. We look for a classical BH solution defined by a spatial discontinuity at the event horizon. This requires non-static solutions with radial fluid velocity ≠ 0 relative to the outside SW observer. The two key questions we want to address here are: What are possible metrics for the inside of such a physical BH? What is the meaning of the BH mass measured by an outside SW observer, like us?

We find here a new solution to these questions, which we call the BHU metric. We will also explore the idea that our Universe corresponds to such BHU solution. As the universe expands tends to Λ which corresponds to a trapped surface Λ = 1/ Λ , just like the event horizon of a BH. Moreover, the density of our universe in that limit is = 3 2 Λ /8 which corresponds to that of a BH, in Eq.1 for = Λ . This indicates that we actually live inside a very massive physical BH. It also tells us what is the metric inside a BH: our Universe is the only object whose interior we know and has the density of a BH. We will show that such situation is a solution to classical GR without singularities. This seems to avoid the singularity theorems in the lines pointed out by Senovilla (1998). But exact mathematical rigour is not the scope or concern of our paper, but rather a more physical and empirical approach to these issues. The idea that the universe might be generated from the inside of a BH is not new and has extensive literature. Easson & Brandenberger (2001) and Oshita & Yokoyama (2018) present a good summary of past and recent literature which mostly focused in deSitter (dS) metric with a dual role of the BH interior and an approximation for our universe. Many of the formation mechanisms involve some modifications or extensions of GR (e.g. Zhang 2018), often motivated by quantum gravity or string theory. This is what we try to avoid here (Ellis & Silk 2014). There are also some examples (e.g. Daghigh et al. 2000) which presented models within the scope of a classical GR and classical field theory with FV interior similar to our BH.fv solution here. These models are affected by the no-go theorem, such us Galtsov & Lemos (2001), that state that no smooth solution to ( ) can interpolate between dS and SW space. But this is not an issue for our solution for two reasons. First, the external asymptotic space is really SW+dS or FLRW (a BH is a perturbation within a FLRW metric), where solutions do exist (e.g. Dymnikova 2003). Second, we do not need ( ) to smoothly transit between metrics: ( ) is trapped in a FV, which is discontinuous by nature.

These previous works provide support to the idea that our universe could indeed be inside a BH, but often they are too simplistic, as they don't contain any matter or radiation. There are some papers that proposed the FLRW metric as BH interior. [START_REF] Pathria | [END_REF] found that for a SW+dS metric behaves like in a closed ( = +1) FLRW metric, but Knutsen (2009) rightly pointed out the inconsistencies in notation and such interpretation. Stuckey (1994) found that a dust dominated FLRW universe could be mathematically joined to an outside BH.SW metric. This is a good precedent to our BH.u solution (in §4.2), but it did not include a Λ term, a FV vacuum, radiation or a physical interpretation. Our ideas are also quite different from Smolin (1992) who speculated that all final (e.g. BH) singularities 'bounce' or tunnel to initial singularities of new universes. Here we propose the opposite, that such mathematical singularities are not needed to explain the physical world. As stated by Ellis (2008), the concept of physical infinity is not a scientific one if science involves testability by either observation or experiment. The BHU model avoids the initial causal and entropy paradoxes (Dyson et al. 2002;[START_REF] Penrose | Conf. Proc. C[END_REF]) because of its origin within a larger spacetime.

In §2 we present our notation for Einstein field equations for perfect fluid and two homogeneous solutions: a FV and an expanding FLRW universe. In §3 we present a brief introduction to the general case of in-homogeneous solutions with spherical symmetry in proper SW coordinates. The FLRW solution can also be expressed in these coordinates. This duality is a key ingredient to find our solution for a physical BH interior. As far as we know this is a new result. In §4 we present two new solutions for the physical BH interior. In §5 we discuss how to apply these solutions to our universe. We end with a summary and a discussion of observational windows to test the BHU.

HOMOGENEOUS SOLUTIONS

We will solve Einstein's field equations (EFE) [START_REF] Padmanabhan | Gravitation[END_REF]:

+ Λ = 8 ≡ - 16 √ - ( √ -L) , (2)
where ≡ -1 2 and L is the matter Lagrangian. For perfect fluid in spherical coordinates:

= ( + ) + (3)
where is the 4-velocity ( = -1), , and are the energymatter density and pressure. This fluid is in general made of several components, each with a different equation of state =

. For a fluid moving with relative radial velocity with = ( 0 , , 0, 0), we have 2 0 = -00 (1 + 11 2 ) and:

0 0 = --2 ( + ) 11 ; 1 1 = + 2 ( + ) 11 1 0 = ( + ) 0 ; 2 2 = 3 3 = (4)
For a comoving observer = 0. The outside manyfold M + is empty space so the solution + is the BH.SW. Because the inside M -is causally disconnected, M + acts like a simple boundary condition (Gaztañaga 2021b). Given some and inside , we will solve EFE inside with such boundary condition to find -. To impose the boundary at we will use the same (proper) SW coordinate frame that is not moving with the fluid so that 1 0 ≠ 0. This could result in solutions for M -that are not static. We will then verify Israel (1967) conditions to double check that the join manyfold M = M -∪ M + is also a solution to EFE and there are no surface terms (see §4.3). This is different from just matching two metrics. We will also find inside and outside the physical BH. This is not enough to find a physical solution to EFE. We also need to study its stability and find a causal formation mechanism.

Scalar field in curved space-time

Consider a minimally coupled scalar field = ( ) with:

L ≡ -= - 1 2 -( ) (5)
The Lagrange equations are: ∇2 = / . We can estimate ( ) from its definition in Eq.2 to find:

( ) = + ( -) (6)
comparing to Eq.3:

= + ; = | | - (7)
In general we can have ≠ ⊥ for non canonical scalar fields (see Eq.5 in Díez-Tejedor & Feinstein 2006 for further details). The stable solution corresponds to = -≡ -:

∇2 = / = 0 ; ≡ = -= ( ) = (8)
where is trapped in the true minimum 0 or some false vacuum (FV) state = 0 + Δ. The situation is illustrated in Fig. 1. The solution to Eq.2 for constant = -= (without matter or radiation) for a general metric with spherical symmetry in proper coordinates (i.e. Eq.15) is given by dS metric in Eq.21 with 2 Λ ≡ 8 Λ /3 where Λ = +Λ/8 . This metric is static which indicates that the vacuum solution is in equilibrium. Quantum tunnelling can result in a phase transtion or vacuum evaporation, which could break the equilibrium or result in slow evolution.

The FLRW metric in comoving spherical coordinates

The FLRW metric in comoving coordinates = ( , , , ), corresponds to an homogeneous and isotropic space-time:

2 = = -2 + ( ) 2 2 + 2 2 (9)
where we have introduced the solid angle: ≡ sinc( √ ) with 2 = cos 2 2 + 2 and is the curvature constant = {+1, 0, -1}. For the flat case ( = 0) we have 2 = 2 . The scale factor, ( ), describes the expansion/contraction as a function of -0 with respect to the true vacuum at 0 = ( 0 ) (blue circle). A FV can be trapped inside other FV. Some FV are unstable or slow rolling, like 7 . Quantum tunneling (dashed green lines) allows to jump between FV. A FV could also be reached by adding kinetical energy to the system to overcome the energy barrier where FV are trapped.

comoving or cosmic time (proper time for a comoving observer). For a comoving observer, the time-radial components are:

00 10 01 11 = ( ) 0 0 ( ) 2 (10)
i.e. = 0 in Eq.4. The solution to EFE in Eq.2 is:

3 = = Λ -4 ( + 3 ) (11) 2 ≡ 2 = 2 0 Ω -3 + Ω -4 + Ω -2 + Ω Λ (12) Λ ≡ vac + Λ 8 (13) ≡ 3 2 8 ; Ω ≡ ( = 1) (14)
where Ω (or ) represent the matter density today ( = 1), Ω is the radiation, vac represents vacuum energy: vac =vac = ( ) in Eq.8, and Λ = -Λ is the effective cosmological constant density. Given and at some time, we can use the above equations to find = ( ) and determine the metric in Eq.9.

Eq.11 also corresponds to the geodesic acceleration (see Eq.43 bellow) and is valid for a perfect fluid in a more general metric. Note how < -/3 (or Λ > 0) produces acceleration (a repulsive gravity) and therefore expansion. This will be the physical mechanism that will avoid the BH singularity for the BHU solution.

During inflation, was dominated by a vacuum field so that Λ ( ) ≡ 32 Λ /8 , which results in Λ . Recent observations show that the expansion rate today is also dominated by Λ . 2 This indicates that the FLRW metric lives inside a trapped surface 1/ Λ = (8 Λ /3) -1/2 , which behaves like the interior of a BH.

2 This is no coincidence: the current late time cosmic acceleration can be understood as caused by inflation, see §5.1.

PROPER COORDINATES

The most general shape for a metric with spherical symmetry in proper or SW coordinates ( , , , ) is Padmanabhan (2010):

1 + 2Φ + 2 2 (15)
where was introduced in Eq.9 to allow for non-flat space. Ψ( , ) and Φ( , ) are the two gravitational potentials. The Weyl potential Φ is the geometric mean of the two:

(1 + 2Φ ) 2 = (1 + 2Φ)(1 + 2Ψ). ( 16
)
Ψ describes propagation of non-relativist particles and Φ the propagation of light. For =we usually get Ψ = Φ = Φ . Eq.15 can also be used to describe the BH.SW solution (or any other solution) as a perturbation (2|Φ| < 1) around a FLRW background:

2

-(1 + 2Ψ) 2 + (1 -2Φ) 2 2 + 2 2 2 (17)
where = ( ) and . This same expression follows from perturbing the FLRW metric in Eq.9.

EFE for the general case of Eq.15 are well known, e.g. see Eq.(7.51) in [START_REF] Padmanabhan | Gravitation[END_REF]. The case =results in 0 0 = 1 1 , and the solution is:

Φ = Ψ = -/ ∫ 0 ( ) 4 2 -Λ 2 /6 -1 / ( 18
)
where 1 is an integration constant. The remaining EFE, 2 2 = 3 3 are equivalent to energy conservation ∇ = 0. For = 0:

= - + 1 + 2Φ Φ. ; = + 1 + 2Ψ Ψ (19)
Note how =results in constant and everywhere, but with a discontinuity at 2Φ = 2Ψ = -1 ( = ). This allows and to be constant, but different in both sides of 2Φ = 2Ψ = -1. One way to address this more formally is to study junction conditions (see §4.3). Another way around this, is to consider anisotropic pressure ≠ ⊥ (Brustein & Medved 2019;Dymnikova 2019) which can result from non canonical scalar field Díez-Tejedor & Feinstein (2006). Empty space ( = = Λ = 0) in Eq.18 results in the BH.SW metric:

2Φ = 2Ψ = -2 / ≡ - / (20)
where 1 = . In the presence of Λ, the BH.SW metric corresponds to Λ = ( ) + Λ/8 = 0 in Eq.8. Outgoing radial null geodesics cannot leave the interior of , while incoming ones can cross inside. The solution to Eq.18 for = = = 0, but Λ ≠ 0 is the deSitter (dS) metric:

2Φ = 2Ψ = -2 2 Λ ≡ -2 / 2 Λ ( 21
)
which is also static and has a trapped surface at = Λ (2Φ = -1) where 2 Λ ≡ 8 Λ /3 and Λ = Λ/(8 ) + ( ). This metric corresponds to the surface of an sphere in a higher dimensional flat spacetime with a constant Ricci curvature = 4Λ. This corresponds to a repulsive force, while Eq.20 is attractive (see Eq.43). The inside of Λ is causally disconnected, like the FLRW metric. Both metrics are equivalent for = Λ (see Mitra 2012) which explains why the dS metric reproduces primordial inflation. The Steady-State Cosmology (SSC), with a perfect cosmological principle, is also reproduced by the same dS metric (O'Raifeartaigh & Mitton 2015). But contrary to the original SSC proponents (Bondi & Gold 1948;Hoyle 1948), Λ = ( ) is trapped to a fixed value and there is no need for continuous creation of matter (C-field). Λ and different values of . All quantities are in units of 1/ Λ . For early time or small :

. A fix acts like an Horizon: as ⇒ ∞ we have ⇒ln (dotted), which freezes inflation to: = ⇒ -ln ( Λ ) = 1/ Λ (dashed).

When we express = with = 4 3 /3 and = Λ , the BH.SW metric in Eq.20 becomes dS in Eq.21:

2Ψ = 2Φ = -2 / = -8 Λ 2 = -2 Λ 2 (22)
We will see this equation again soon (e.g. Eq.39). In comoving coordinates, dS singularity corresponds to a comoving Hubble horizon that shrinks to zero (see Fig. 4). But note that this singularity can not be reached from the inside. Radial null events ( 2 = 0) connecting (0, 0 ) with ( , ) follow:

= + 0 -0 2 / -1 / + 0 -0 2 / + 1 (23)
so that it takes = ∞ to reach = from any point inside. So the inside observer is trapped as in the FLRW case.

When we have both and Λ constant everywhere, the solution to Eq.18 is: 2Φ = 2Ψ = -2 2 Λ -/ , which corresponds to SW-dS metric, a BH in a dS background which has slightly shifted event horizon and is an approximation to our FLRW metric today. where Λ ≡ is the Jacobian. Using the duality transformation :

The FLRW metric in proper coordinates

Λ = = (1 + 2Φ ) -1 (1 + 2Φ ) -1 , (24)
it is straightforward to verify that Eq.15, with 2Φ = -2 2 and arbitrary ( ) and Ψ, transform into the FLRW of Eq.9:

= Λ -(1 + 2Ψ) 0 0 (1 + 2Φ) -1 Λ = -1 0 0 2 , ( 25
)
In other words, these two metrics are the same:

-(1+2Ψ) 2 + 2 1 -2 2 + 2 2 = -2 + 2 2 + 2 2 (26)
So a generic spherically symmetric in-homogeneous metric can look like the homogeneous FLRW metric in comoving coordinates. And a general FLRW metric can be cast in proper coordinates. The expression for Ψ = Ψ( , ) and = ( , ) or = ( , ), depends on ( ). For ( ) = Λ (i.e. = Λ ) we have Ψ = Φ and

= ( , ) = - 1 2 Λ ln [1 -2 Λ 2 2 ], (27)
where < Λ = 1/ Λ , which reproduces dS metric. In comoving coordinates the metric is inflating exponentially: = Λ , while in proper coordinates it is static. Fig. 2 illustrates how this is possible and shows = ( , ) (see Mitra 2012 for some additional discussion). The remarkable frame duality of Eq.26, from a comoving time frame to a proper SW frame, is a new result as far as we know, and a key ingredient to interpret our new physical BH solution.

BLACK HOLE SOLUTIONS

False Vacuum Black Hole (BH.fv) solution

Eq.20 and Eq.21 are the simplest solutions to EFE. They both correspond to some form of empty space. We can use the same physical models to explore the simplest physical BH interior:

= -= 0 for > 0 + Δ for < (28)
where Δ > 0. To recover BH.SW solution outside, we need Λ = -8 0 so that Λ = 0. But in a realistic situation, on larger scales the BH.SW metric should be considered a perturbation of FLRW background, e.g. Eq.17. It is straightforward (but tedious) to present the solution in this way. We take 1 = 0 in Eq.18 to avoid singular solutions. The solution to Eq.18 (which we called BH.fv) is then:

2Φ = 2Ψ = - / for > ≡ 2 -2 2 Λ for < = Λ ≡ 1/ Λ (29)
where: Λ = = Δ and = 4 3 3 Δ. This has no singularity at = 0. Note how, contrary to what happens in the BH.SW, in the BH.fv solution, the metric components don't change signature as we cross inside . In both sides of we have constant but different values of and . This comes from energy conservation in Eq.19. There is a discontinuity ( = ∞ and = ∞) at 2Φ = -1 where = , in agreement with Eq.19, but the metric is static and continuous at . This solution only happens when = (8 Δ/3) -1/2 . The smaller Δ the larger and more massive the BH. In the limit Δ ⇒ 0, we have ⇒ ∞ and we recover Minkowski space, as expected.

At a fixed location, the scalar field inside the BH is trapped in a stable configuration ( = 0 + Δ) and can not evolve ( = 0 in Eq.7). The same happens for the field outside (see Fig. 1). A FV in Eq.28 with equal Δ but with smaller initial radius = < is subject to a pressure discontinuity at = which is not balanced in Eq.19 and results in a bubble growth (Aguirre & Johnson 2005). Such boundary grows and asymptotically reaches = (see top panel of Fig. 2, Fig. 3 and Eq.41). The inside of is causally disconnected, so the pressure discontinuity does not act on = , which corresponds to a trapped surface.

Black Hole Universe (BH.u) solution

We will look next for solutions to the metric of Eq.15 where we also have matter = ( , ) and radiation = ( , ) inside:

( , ) = -= 0 for > 0 + Δ + + for < (30) with = -0 -Δ + /3 ≠ -inside.
This means Φ ≠ 0 and ≠ 0: the fluid inside has to move relative to SW frame of the outside observer3 . For > , the solution is the same as Eq.29. For the interior we define: 2Φ ≡ -2 2 ( , ), so that:

2Φ( , ) = - / for > = 2 -2 2 ( , ) for < = 1/ Λ (31)
where = 2 = 1/ Λ as before. We can find the interior solution with a change of variables of Eq.24-26. The solution to Eq.31 is then ( , ) = ( ). Given some ( ) and ( ) in comoving coordinates in the interior of a BH we can use Eq.12 with Λ = Δ = 3 -2 /8 to find ( ) and ( ). We call this a BH universe (BH.u). To complete the solution, i.e. to find Ψ and = ( , ), we need to solve Eq.24 with 2Φ = -2 2 ( ). For ( ) = Λ the solution is Ψ = Φ and Eq.27. The flat FLRW metric with = Λ becomes dS metric in Eq.21 as in the BH.fv solution.

Given in Eq.10 we can find ¯ in the proper frame using the inverse matrix of Eq.24: ¯ = (Λ -1 ) (Λ -1 ) :

¯ 0 0 = - -2Φ 1 + 2Φ ; ¯ 1 1 = -2Φ 1 + 2Φ (32)
which is independent of Ψ. Comparing to Eq.4 gives the velocity in the proper frame 2 = -2Φ = 2 2 , which is just the Hubble law. The Lorentz factor is = (1 + 2Φ) -1/2 so that gives the proper length, in agreement with Eq.15.

Solution ( , ) = ( ) in Eq.31 is valid for all < = 1/ Λ because ( ) > Λ . We can see this by considering outgoing radial null geodesic in the FLRW metric of Eq.9:

= ( ) ∫ ∞ ( ) = ∫ ∞ 2 ( ) < 1 Λ = (33)
which shows that signals can not escape from the inside to the outside of the BH.u. But incoming radial null geodesics ( ) ∫ 0 ( ) can in fact be larger than if we look back in time. This shows that inside observers are trapped inside the BH.u but they can nevertheless observe what happened outside (Gaztañaga & Fosalba 2021).

Junction conditions

We can arrive at the same BHU (BH.fv and BH.u) solutions using the junction conditions of Israel (1967). Here we follow closely the notation in §12.5 of [START_REF] Padmanabhan | Gravitation[END_REF]. We will combine two solutions to EFE with different energy content, as in Eq.30, on two sides of a timelike hypersurface Σ = M -∩ M + . The inside -is FLRW metric (or dS metric for = Λ ) and the outside + is BH.SW metric. This is similar to the case §12.5.1 in [START_REF] Padmanabhan | Gravitation[END_REF] with the difference that we use = 0 (instead of = 1) and consider a general FLRW solution ( ) with Λ, and (instead of a pressure-free dust model without Λ.). This is relevant to provide the limiting trapped surface = 3/Λ. We define Σ to be fixed in comoving coordinates at = , so Σ only depends (here we fix = 1 when =

). For the outside SW coordinate system, Σ + is described by = ( ) and = ( ), where is the comoving time in the FLRW metric. We then have:

= ; = , (34)
where the dot refers to derivatives with respect to . The induced metric ℎ -on the inside of Σ -with = ( , , ) and fixed = , is:

2

Σ -= ℎ - = -2 + 2 ( ) 2 2 ( 35
)
has to agree with ℎ + , the BH.SW metric outside at Σ + :

- 2 + -1 2 + 2 2 = -( 2 -2 / ) 2 + 2 2 (36)
where = 1 -. The matching condition ℎ -= ℎ + is:

( ) = ( ) ; = 2 + ≡ ( , ) (37)
Thus, for a given FLRW solution ( ) we know both and . The extrinsic curvature ± normal to Σ ± from each side is:

-= 0 ; -= -= - 1 + = ; + = + = - (38)
Thus, the second matching condition -= + requires = 1, which using Eq.37 results in:

2 = 2 2 = (39) = 1 1 -2 2 ⇒ = ∫ (1 -2 2 ) (40)
This results in 2Ψ = 2Φ = -2 2 = -/ in the junction Σ (same as Eq22 for = Λ ). This is a generalization of dS space for arbitrary ( ): 2Ψ = 2Φ = -2 / 2 with ≡ 1/ ( ). The critical density inside , corresponds to that of a BH: from Eq.39 we have 2 = / 3 = 8 /3. This Σ junction grows and tends to = 1/ Λ . It takes = ∞ in the SW time of Eq.40 to asymptotically reach = 1/ Λ (see Fig. 2). In this limit, Eq.39 reproduces the BH.u junction of Eq.31 (and the BH.fv junction of Eq.29 for constant = Λ ). Before that, the BHU junction is not static (not even in the SW frame) as decays into Λ when grows to . Despite the discontinuity in at , the BHU metric and extrinsic curvature are continuous when we join them with the expanding timelike hypersurface of Σ. This proofs that the BHU metric is also a solution to EFE and there are no surface terms in the junction (see Eq.21.167 in [START_REF] Misner | Gravitation Mitra A[END_REF]). This does not require that both of the joined metrics have identical Riemann tensor or associated invariant scalars. 4 The Λ term corresponds to a trapped surface = 1/ Λ in the FLRW (or dS) metric which matches the horizon of a BH in empty space (see Fig. 5).

Evolving junction: internal BH dynamics

The junction conditions indicate that the division between interior and exterior solutions in Eq.29 and Eq.31 is not necessarily , which is the limiting case. This is illustrated in Fig. 3. That both the metric and the external curvature are continuous at shows that there are no surface terms and the join metric is a solution to EFE (see Eq.21.167 in [START_REF] Misner | Gravitation Mitra A[END_REF]). The energy-momentum tensor corresponding to this solution has a discontinuity (as expected for a BH): Λ = -Λ ≠ 0 for < and Λ = 0 for > .

Inside the physical BH we have an expanding junction: = ( ) = ( )

. Because ( ) < is always inside , the external SW observer can not distinguish this evolving junction from the limiting static one =

. This is why we chose to express the solution this way. The junction ( ) grows and asymptotically tends to as shown in Fig. 3. This happens at a finite comoving time Λ as in the top panel of Fig. 2. The exact function depends on the form of = ( ). For constant / = = Λ = 1/ , the solution can be expressed analytically as:

( ) = 0 = Λ ( -Λ ) (41)
where we have chosen = 1 when = . We start with a finite size = 0 = 0 at = 0, where 0 = -Λ Λ . After Λ Λ e-folds, 0 grows into = . This inflation stops asymptotically at = Λ =ln 0 . We can think of 0 as a quantum size (FV) fluctuation of energy Λ = Δ, which (in empty space) will inflate to size = 1/ Λ = (8 Δ/3) -1/2 . This new solution to EFE is not just an arbitrary matching of two other random solutions. It is a new solution of a new physical configuration given by the energy content in Eq.30. This configuration corresponds exactly to our definition of a generic physical BH. The one we set to find in the introduction and whose horizon separates two regions with different matter-energy content. The same horizon defines the junction of two well known solutions to EFE.

IMPLICATIONS FOR OUR UNIVERSE

The BH.fv interior, dS metric, can be transformed into a FLRW metric with constant = Λ . This frame duality provides a new interpretation for the BH.fv solution in Eq.29. This is not only a solution for a BH inside a universe. The inside comoving observer, sees this solution as an expanding inflationary universe inside a BH, even when the metric is static in proper coordinates and = is fixed. The same happens with the BH.u solution of Eq.31, which is equivalent to a child FLRW in the interior.

Recall how the outside BH.SW solution should be considered a perturbation of a parent FLRW in Eq.17. So we have two nested FLRW metrics which are connected with the BHU. Each one could have a different effective Λ (or FV). So could our universe be a child FLRW metric? The fact that we have measured Λ ≠ 0 provides a strong indication that this is the case. It is hard to explain what Λ means or the coincidence problem otherwise (see below).

The change of variables in Eq.24 is only valid for proper coordinates that are centered at the center of the BH location. But in the transformed (comoving) frame of Eq.26 any point inside is subject to the same expansion law with equal ( ). From every point inside de BHU, comoving observers will see an homogeneous and isotropic space-time around them. At least all the way back in time to = 1/ Λ , which is in their distant past. All points inside the BH were at the center in their distance past (just as in the homogeneous expanding universe). Even if all points inside see an homogeneous expansion, as they look back in time from a position that is off centered, some regions of the sky will be closer to the past horizon than others. This could result in some significant deviations from isotropy and homogenety on the largest scales. Regions outside the trapped surface could see a similar background but with uncorrelated fluctuations that fit different parameters. Such deviations have already been measured in the form of large scale CMB anomalies and variations of cosmological parameters (see Fosalba & Gaztañaga 2021).

Note how we can have FVs inside other FVs (see Fig. 1). So we can have BHs inside other BHs or FLRW metrics inside other FLRW universes. Mathematically this looks like a Matryoshka (or nesting) doll or a fractal structure. But physically, in the common SW frame, each BH has a different mass and therefore different physical properties. The child FLRW BHU have smaller mass (and larger FV) than the parent BHU. A BHU of one solar mass can have a FLRW metric inside but this inside will not have any galaxies and is going to be very different from that in a 5.8 × 10 22 BHU, like ours, which contains billions of galaxies and BHs of many different sizes. So each BHU layer could be physically quite different from the next, unlike Matryoshka dolls or fractal structures.

The evolution of the BH universe

How did the universe evolve into the solution of Eq.31? This is an important question. It is not enough to find a solution to EFE. We need to make sure that such a configuration can be achieved in a causal way. Without Λ, the FLRW universe has no causal origin: the Hubble rate (in Eq.12) is the same everywhere, not matter how far, and this is not causally possible. The comoving coordinate = that fixes the junction in §4.3 above can be identified as the causal horizon § in the zero action principle (Gaztañaga 2021b). In the FLRW Universe, the Hubble Horizon is defined as = / . Scales larger than cannot evolve because the time a perturbation takes to travel that distance is larger than the expansion time. This means that > scales are "frozen out" (structure can not evolve) and are causally disconnected from the rest. Thus, / represents a dynamical causal horizon that is evolving.

The standard evolution of our universe is shown in Fig. 4. Note that here we choose = 1 now, as opposed to Fig. 3 where = 1 corresponds to =

. It turns out that both are not so different (the so call coincidence problem in cosmology). A primordial field settles or fluctuates into a false (or slow rolling) vacuum which will create a BH.fv with a junction Σ in Eq.35, where the causal boundary is fixed in comoving coordinates and corresponds to the particle horizon during inflation § = /( ) or the Hubble horizon when inflation begins. The size = ( ) § of this vacuum grows and asymptotically tends to = / following Eq.39 with =

. The inside of this BH will be expanding exponentially = while the Hubble horizon is fixed 1/ . Accooring to standard models of primordial inflation (Starobinski ǐ 1979;Guth 1981;Linde 1982;Albrecht & Steinhardt 1982), this inflation ends and vacuum energy excess converts into matter and radiation (reheating). This results in BH.u, where the infinitesimal Hubble horizon starts to grow following the standard BB evolution.

Note that the inflation in the BH.fv solution (i.e. Eq.41) stops naturally at cosmic time = --1 ln § (see Fig. 2) when proper distance is = ( ) § = 1/ . In standard models of primordial inflation, is much larger that Λ so that 1/ is much smaller than 1/ Λ . So a quantum FV fluctuation Δ only grows to a maximun size = = (8 Δ/2) -1/2 = 1/ . Something else has to happen if we want the size to become cosmological. It could be reheating or some other mechanism. Quantum tunneling into smaller Δ (see Fig. 1) also produces larger . Matter and radiation can also appear some other ways: from the original quantum fluctuations, quantum tunneling from/to other FV or infall of matter from outside (see §6.3-6.4). Regardless of these details, § remains the causal scale for the original BH.fv inflation in Eq.41. Recall that the BH.fv solution requires a discontinuity in Λ = 0, so this BH.fv evolution happens with independence of what we assume about Λ in EFE. A causal boundary in empty space generates a boundary term in the action that fixes the value of Λ to Λ = 4 < + 3 >, where the average is over the light-cone inside § (Gaztañaga 2021b). This Λ represents a trapped surface for the emerging BH.u universe. 5The observable universe (or particle horizon) after inflation is:

= ( ) = ∫ ln ( ) = (1) -¯ ( ), (42)
where is the scale factor when inflation ends. For Ω Λ 0.7, the particle horizon today is

(1) 3.26 / 0 and ¯ ( ) = ∫ 1 ln /( ) is the radial lookback time, which for a flat universe agrees with the comoving angular diameter distance, = ¯ . The observable universe becomes larger than = Λ when > 1, as shown in Fig. 4 (compare dotted and dashed lines). This shows that, observers like us, living in the interior of the BH universe, are trapped inside but can nevertheless observe what happened outside. We can estimate § from Λ =< /2 + >, where the average is in the lightcone inside § . For Ω Λ 0.7 Gaztañaga (2021b) found: § 3.34 / 0 which is close to today. But imagine that Ω Λ is caused by some DE component and has nothing to do with § . We still have that § , because otherwise § would have cross = 1 early on, resulting in smaller than measured (see Fig. 4). Thus, at the time of CMB last scattering (when ), § corresponds to an angle = § / 1 rad 60 deg. So we can actually observe scales larger than § . Scales that are not causally connected! This could be related to the so-called CMB anomalies (i.e, apparent deviations with respect to simple predictions from ΛCDM, see Planck Collaboration 2020b and references therein), or the apparent tensions in measurements from vastly different cosmic scales or times (e.g. Planck Collaboration 2020a).

DISCUSION & CONCLUSION

We have looked for classical non-singular GR solutions for the BH interior. Our motivation is find a physical model and study if this results in some different properties for observed BHs. The outside manyfold M + of a BH is approximated as empty space so the solution + is the BH.SW metric. Because the inside M -is causally disconnected, M + acts like a simple boundary condition. Given some and inside , we have solve EFE inside with such boundary condition to find -, the inside metric of a physical BH. To our surprise we have found that -is just the well known FLRW, the same metric that describes our universe! This frame duality, represented by Eq.24, has several observational consequences, as we will discuss below.

To impose the boundary at we have use the same (proper) SW coordinate frame that is not moving with the fluid so that 1 0 ≠ 0. This results in a solution for M -that is not static. We have verified Israel (1967) conditions to double check that the join manyfold M -∪ M + is also a solution to EFE and there are no surface terms (see §4.3). This is different from just matching two random metrics. We have find inside and outside the physical BH and we have found what the BH mass means (see §6.5 below).

A key point in these new solutions is the concept of FV energy (see §2.1 and Fig. 1) and the fact that there is a FV or Λ discontinuity at = (which could also be understood as a boundary to the Einstein-Hilbert action, Gaztañaga 2021b). The relativistic version of Poisson equation comes from the geodesic deviation equation (see Eq.12 in Gaztañaga 2021b):

∇ g = Θ + 1 3 Θ 2 = = Λ -4 ( + 3 ) ( 43
)
where g is the geodesic acceleration [START_REF] Padmanabhan | Gravitation[END_REF]. This is the same as the Raychaudhuri equation for a shear free, non-rotating fluid where Θ = ∇ and is proper time. The above equation is purely geometric: it describes the evolution in proper time of the dilatation coefficient Θ of a bundle of nearby geodesics. This is valid for a perfect fluid in a general metric. Note how < -/3 (or Λ > 0) produces acceleration (a repulsive gravity) and therefore expansion. This is the physical mechanism that explains how a FV (which has = -) avoids the classical singularity theorems in GR (see Senovilla 1998 and references therein). We can add both matter and radiation to both sides of and we still have a BHU solution. The BHU connects two FLRW metrics (see Fig. 5).

False Vaccum BH solution (BH.fv)

BH.fv corresponds to constant FV discontinuity (Eq.28) with dS metric inside (Eq.29). dS metric has a trapped surface at = which matches the BH.SW event horizon. A constant density (or negative pressure) corresponds to a centrifugal force, 2Φ = -( / ) 2 that opposes gravity, 2Φ = -/ , i.e. Eq.18. The equilibrium happens when both forces are equal, which fixes = , and correspond to stable circular Kepler orbit. The BH mass is just = Δ. Similar solutions have been known for a while (e.g. Easson & Brandenberger 2001;Daghigh et al. 2000;Firouzjahi 2016;Oshita & Yokoyama 2018;Dymnikova 2019), such as the gravastar (Mazur & Mottola 2015) or bubble/baby universes (Aguirre & Johnson 2005;Kusenko & et al. 2020). The BH.fv solutions here is not just a mere matching of the BH.SW and dS metrics, but a proper solution to Einstein's Field Equations (EFE) that fits our definition of a physical BH. We have also shown in §4.3 that there are no surface terms in the junction, Σ, which is a timelike expanding hypersurface. This is different from the gravastar model (which has matter in the surface) or the anisotropic models with negative radial pressure (Brustein & Medved 2019;Dymnikova 2019).

The BH universe solution (BH.u)

In our second solution (Eq.31), the BH interior is the FLRW metric. This BH.u solution is new, as far as we know. As discuss in the introduction, previous attempts were not proper solutions [START_REF] Pathria | [END_REF]Knutsen 2009) or did not include radiation or Λ (Stuckey 1994). We can have other BHs, matter and radiation inside a BHU within a larger space-time. The inside needs to be expanding as in the FLRW metric of Eq.9, with the same trapped surface given by Λ . This holds the expansion and balance gravity at as in the BH.fv solution. The join solution (Eq.31) is also a solution to Einstein's field equations as the two metrics reduce to the same form on a junction of constant = in Eq.35, and the extrinsic curvature in Eq.38 is the same in both sides. The junction conditions indicate that the division ( ) between interior and exterior solutions in Eq.29 and Eq.31 is not necessarily , which is the limiting case. The junction ( ) grows and asymptotically tends to as shown in Fig. 3. The exterior metric could also be FLRW, as the BH.SW metric can be considered a local perturbation within a larger FLRW background with arbitrary = ( ) in Eq.17. This is illustrated in Fig. 5. In this case we need to distinguish between two different effective Λ , the one in the inside FLRW metric, Λ = Δ and the one for the outside background, which should be smaller.

The solutions to the field equations are independent of the choice of coordinates but ¯ ( , ) depends on the fluid motion (see Eq.32). We used comoving coordinates ( , ), where the fluid is expanding and the observed is comoving, to find the interior solution. But we can then transform back to proper SW frame ( , ), using the duality transformation Eq.24, to find a full BH solution in Eq.31 that is continuous in the metric and curvature at , like in the BH.fv case. As in the singular BH.SW metric, outgoing radial null geodesics cannot escape the event horizon, but incoming ones can enter (see discussion around Eq.33). So the BHU solution is a physical BH.

BH formation

Another issue, which we only address partially here, is how such physical BH solutions can be achieved (e.g. astrophysical and primordial BH formation) and if they can have a causal origin. A local (quantum) fluctuation Δ could generate a physical BH solution because either jumps into a FV or because of the zero action principle (see the Appollonian Universe in Gaztañaga 2021b). Hubble dumping of the kinetic energy of a classical scalar field (see Fig. 1) can also result in a FV trapped field configuration. Such initially small local discontinuity, with FV energy density Δ, will grow as Eq.41 until it reaches BH size corresponding to = Δ. So Δ is the BH density: the smaller Δ the larger the BH size and mass. Can we then make a very massive BH out of a quantum fluctuation, for free? There are some caveats to this. First, quantum fluctuation have a size given by ℏ so typical values of Δ are much larger than Λ in our universe. Second, if Δ is small, there is a chance that could tunnel back to the true vaccum (see Fig. 1) before the BH grows to it's SW radius .

Another caveat is that we have obtained these physical BH solutions for empty space outside (i.e. when Φ is a small perturbation within a FLRW background of Eq.17). When and become large compared to the background these solutions are not valid.

Could a BHU formed just from the final collapse of a dying star? Instead of a forming a BH singularity (as usually assumed), such collapse may just result into a large supernova (SN) explosion, where most of the stellar mass disappears in the SN and only a white draft or neutron star remains. Also a small BH could remain but it can not be a BH.SW because this is non physical. If such small BH is a BHU, quantum tunneling into smaller values of Δ will make the BH grow as ∝ Δ -1/2 . The SN explosion could also provide some kinetical energy = Δ needed for to jump into a higher FV (see Fig. 1), which will reduce . This will trigger a BH.fv formation (so a new expansion in Eq.41). All these are speculation as we need more detail modeling to find how this happens. But the point we want to make is that the BH interior is important for models of BH formation and we can not just assume that a BH is a singular BH.SW metric inside, because this is not physical. In Fig. 1 we can see that BHU of very different masses can be formed as the spectrum of Δ values could be quite broad if we allow ( ) to be a superposition of many plane waves. But, as we will show next, it is unlikely that BH mass can grow from accretion, as it is usually assumed for BH.SW.

Gravitational Waves & Cosmic rays

Matter and radiation can infall inside the BHU solution. This results in a jump to the internal Hubble expansion rate ( Δ / = 1 /(2 2 )). Such jump will be diluted away by the same internal expansion. This process is anisotropic by nature. So whatever we feed a BHU with, it is converted into internal kinetic expanding energy. For the outside SW frame this loss in the internal energy density, , corresponds to BH mass loss, , which could results in a Gravitational Wave (GW) and high energy cosmic rays, analogous to what happens during the Hubble damping in re-heating. The size of the signal could depend on the BHs spin and the anisotropy of the merging. If a BH or compact object of mass 1 is accreted to a larger BHU, with 2 > 1 , the mass 1 behaves like dark matter. Within a Hubble time, 1 will be diluted away by the internal expansion of 2 , without change to its outside frame observed mass 2 . In a more general case, we could speculate that some of the energy from 1 could be transfer as kinetical energy to the 2 field trapped in 2 , this could also be Hubble damped, as in reheating, into matter and radiation which will then also be diluted by the internal BH expansion of 2 . So again, no changes to 2 . Only when is large enough to overcome the FV barrier (see Fig. 1) the merger could displace to different FV resulting in a new BH mass larger than 2 but likely much smaller than 2 + 1 : 2 2 + 1 . When 1 = 2 with identical FV ( 1 =

2 ) the final mass could be closer to 1 + 2 .

Such BH binaries could provide observational evidence or rule out the BHU solution. In the SW frame, the BH mass 2 will first increase to 2 + 1 and then decrease back to 2 2 + 1 . The upper bound corresponds to a typical mass loss rate in the merger (ignoring GW energy lost from spinnig and ringdown). The mass dilutes as -3 -3 2 , for dS expansion (slower rate when both BH have similar mass), so that:

3 2 -3 2 1 ( 44
)
where is comoving time and -1 2 = = 2 2 is the internal Hubble rate of 2 . The external observer uses the SW frame ( , ).

Using the inverse matrix of Eq.24 and Eq.27 we have:

≡ = - 1 + 2Φ 3 2 -3 2 1 (45)
where we assumed that does not vary with and we are far from . This mass loss has to be compensated with some energy emission in the form of high energy cosmic rays or GW, ( -) . The typical amplitude |ℎ|:

|ℎ| = |ℎ| 2 10 -19 1 Mpc 3 2 -3 2 (46)
where is the luminosity distance. The observed frequency , if emitted at redshift from us, is = 2 (1+ ) so that:

2 ≡ 3 2 2 (1 + ) = 3(2 2 ) -1 2 (1 + ) 100 KHz 2 1 1 + (47)
The characteristic amplitude spectral density strain is then:

2 1/2 |ℎ( )| 10 -19 -1/2 (1 + ) 1 Mpc 2 2 2 2 + 2 1/2 (48)
which is within LIGO sensitivity (Buikema et al. 2020). So this predicts extra power with an exponential decay in the GW signal right after each detected GW ringdown (Carcasona & Gaztanaga, in preparation). This background is stronger at LISA mHz frequencies ( 2 10 7

). This GW background could also be observed in the CMB polarization (from tensor modes).

What is for a physical BH?

The BH mass in the BHU is given by the FV excess energy Δ, so that in Eq.1 is = Δ and = (32 3 Δ/3) -1/2 . So the larger Δ the smaller the BH mass and size. This is independent of the matter and energy content that falls inside the BH (see §6.4). So in the BHU solution does not correspond to the actual total mass or radiation inside, which is not observable from the outside, but should instead be interpreted in terms of the FV energy excess Δ. This could have implications for models of astrophysical BH formation (such as Kormendy & Ho 2013) and primordial BH formation (e.g. Kusenko & et al. 2020 and references therein) which usually assume that BH accretion and merging results in linear increase of the BH mass .

Our universe as a BH

Both BH type solutions can be interpreted as a BH within our universe or as an expanding universe inside a larger space-time. As pointed out in the introduction, that the universe might be generated from the inside of a BH has a long and interesting history. Knutsen (2009) argued that and in the homogeneous FLRW solution are only a function of time (in comoving coordinates) and can not change at = to become zero in the exterior. This is an important point and seems to contradict the BHU solution. The riddle is resolved with Λ. Without Λ the FLRW universe can not have a causal origin: the comoving density and Hubble rate are the same everywhere, and this is not causally possible. A causal horizon § fixes Λ (Gaztañaga 2021b) which solves this problem and also generates an even horizon § = similar to that of a BH.SW: . This allows for an homogeneous FLRW solution inside that has a Λ discontinuity at and looks in-homogeneous in the SW frame. Homogeneity is therefore the illusion of the comoving observer inside Λ =

. The FLRW metric is trapped inside < Λ , and is then equivalent to an inhomogeneous spherically symmetric metric of Eq.26. The FLRW metric is only homogeneous in space, but not in space-time. A new frame where comoving time and space are mixed, can break or restore this symmetry. The frame duality in Eq.24 is only valid for proper coordinates that are centered at the BH location. But in the transformed (comoving) frame any point inside the BHU is subject to the same expansion law with equal ( ). From every point inside de BHU, observers will see an homogeneous and isotropic space-time around them. Just like in the universe around us.

Evidence for a BHU

We can sketch the evolution of our universe with this BHU model (see Fig. 34). In proper coordinates this solution has no BB (or bounce): it is not singular at = 0 or at = 0, because we have a non-singular BH.fv before we start the FLRW BH.u phase. The inside comoving observer is trapped inside < = 2 = 1/ Λ and has the illusion of a BB. The space-time outside (the parent FLRW universe) could be longer and larger than the BB estimates. We could have a network of island universes with matter and radiation in between.

This also explains why our universe (or other island universes) is expanding and not contracting. The initial fluctuation 2 = 8 Δ/3 could be expanding ( > 0) or contracting ( < 0). In the later case it will either recollapse very quickly or it will else bounce into expansion dominated by the repulsive gravitational force that results from the negative pressure from constant Δ or Λ (see Eq.43).

We have other observational evidence that the expanding metric around us is inside a BHU. We can recover the BB homogeneous solution in the limit Δ ⇒ 0, where we have ⇒ ∞ and Λ = 0. But we have measured Λ > 0 (Ω Λ 0.7) which implies 5.8 × 10 22 and / 0 , as in the BHU. The causal interpretation for § , also explains the observed coincidence between Λ and today (Gaztañaga 2020;Gaztañaga 2021b).

If we look back to the CMB times, § corresponds to 60 degrees in the sky. The observed anomalies in the CMB temperature maps at larger scales (Gaztañaga 2020;Gaztañaga 2021b;Fosalba & Gaztañaga 2021;Gaztañaga & Fosalba 2021) provide additional support for the annisotropies expected in the BHU model. There is also a window to see outside our BHU using the largest angular scales for > 2 and measurements of cosmological parameters from very different cosmic times. There is already mounting evidence for this (e.g. Planck Collaboration 2020a; Riess 2019; Abbott et al. 2019).

If there are other island universes outside ours, Galaxies and QSO, as well as BHs, could be accreted from outside Λ into our BHU. Because the horizon 1/ Λ is so large, we can only see evidence of those mergers at early times, during or right after the CMB, when § subtends 60deg. on the sky. Could this be related to rarely old QSO or galaxies observed at high z? If our BHU merges with another BHU which is few % smaller, we might be able to see such % glitches in ( ) with current or future data, at > 2 and very large angle separation. Another possible observational evidence for the BHU solution is outgoing high energy cosmic rays or GW background signal from BH accretion or mergers within our BHU (see §6.4). High energy cosmic rays have been linked with X-ray binaries and AGNs, both hosting BH of different masses. Such GW background signal could also be observable as CMB tensor fluctuations. Camacho & Gaztañaga (2021) found evidence for homogeneity and lack of correlations in the CMB at > Λ . This suggests that the underlying physical mechanism sourcing the observed anisotropy encompasses scales beyond our causal universe. Fosalba & Gaztañaga (2021) found variations in cosmological parameters over large CMB regions. This is the largest reported evidence for a violation of the Cosmological principle. Such observations indicate a breakdown of the standard BB picture in favor of the BHU. Their Fig. 31 shows that the size of these regions follow the BHU relation between § and Λ . This is consistent with the idea that our universe was accreted to or created by a larger BHU.

Figure 1 .

 1 Figure 1. A potential ( ), of a classical field ( ), made of the superposition of 3 plane waves. A configuration with total energy: = + (black dot at 5 ) can slowly loose its kinetic energy 5 (e.g. via Hubble damping) and relax into a static ( = 0) ground state = ≡ ( ), = 1, 2, 3 (red dots), which we call false vacuum (FV) because they have an energy excess Δ ≡-0 with respect to the true vacuum at 0 = ( 0 ) (blue circle). A FV can be trapped inside other FV. Some FV are unstable or slow rolling, like 7 . Quantum tunneling (dashed green lines) allows to jump between FV. A FV could also be reached by adding kinetical energy to the system to overcome the energy barrier where FV are trapped.

Figure 2 .

 2 Figure 2. Logarithm of proper radius = ( ) (top) and comoving time (bottom) as a function of SW time in Eq.27 for ( ) = Λ and different values of . All quantities are in units of 1/ Λ . For early time or small :. A fix acts like an Horizon: as ⇒ ∞ we have ⇒ln (dotted), which freezes inflation to: = ⇒ -ln ( Λ ) = 1/ Λ (dashed).

  Consider a change of variables from = [ , ] to comoving coordinates = [ , ], where = ( ) and angular variables ( , ) remain the same. The metric transforms to = Λ Λ ,

Figure 3 .

 3 Figure 3. Illustration of the interior dynamics of a BHU. The junction = ( ) (black disk) grows towards the SW radius = (dashed red circle). The inside of = ( ) is dS or a FLRW metric (dominated by negative pressure: = -) while the outside is empty ( Λ = 0) BH.SW metric, a perturbation around a FLRW background with other BHs and matter.

Figure 4 .

 4 Figure 4. Comoving (top) and proper (bottom) = ( ) coordinates in units of / 0 as a function of cosmic time (scale factor). The Hubble horizon / (blue continuous line), is compared to the observable universe after inflation (dashed line) and the primordial causal boundary § = /( ) (dot-dashed red line). Larger scales (green shading) are causally disconnected, smaller scales (yellow shading) are dynamically frozen and shrinks to zero in comoving coordinates during dS phases. After inflation / grows again. At 1 the Hubble horizon reaches our event horizon = / Λ .

Figure 5 .

 5 Figure 5. Spatial representation ( = 0) of 2 = (1 + 2Φ) -1 2 + 2 2 for Minkowski (Φ = 0, flat), Schwarzschild (SW, 2Φ = -1/ ), and deSitter (dS, 2Φ = -2 ) metrics. Dashed lines show geodesics. The BHU solution has two nested FLRW metrics and no singularities.

MNRAS 000, 1-10 (2020)

= = -(1 + 2Ψ) 2 + 2

This makes physical sense. Specially because the outside observer has no way to see inside. In general, could also be anisotropic.

This was questioned by two anonymous referees from Physical Review Letters (PRL), who recently rejected a preprint version of this same analysis(Gaztañaga 2021a). Ref-A found that these results warrant publication in PRL. Ref-B said: "I do not think the solutions presented are valid, and I still contend that a valid solution would require peculiar matter at the junction, which would be contrived". Ref-C agreed and added: "the Kretschmann scalar could reveal discontinuities" so the analysis "lacks a sufficiently firm foundation to be useful to research in this area". The solution to Eq.38 in Eq.39 demonstrates that there are no surface terms in Σ. The Kretschmann scalar is different for the BH.SW and FLRW metrics because is different, but this has nothing to do with the junction condition for Σ. ArXiv refused to publish this paper (as well as the preprint version of Gaztañaga 2020; Gaztañaga 2021b) because it "does not contain sufficient original or substantive scholarly research". This of course is just an unjustified opinion (see darkcosmos.com).MNRAS 000, 1-10 (2020)

This picture changes if the fluctuation BH.fv happens within a non empty background. In such case we will just create a BH inside the larger universe, which is dominated by its own Λ and matter contend. See §6.3 below.MNRAS 000, 1-10 (2020)

ACKNOWLEDGMENTS

This work has been supported by spanish MINECO grants PGC2018-102021-B-100 and EU grants LACEGAL 734374 and EWC 776247 with ERDF funds. IEEC is funded by the CERCA program of the Generalitat de Catalunya.

DATA AVAILABILITY STATEMENT

The data and codes used in this article will be shared on request.


Corentin Bernard

Jocelyn Monnoyer

Sølvi Ystad

Michaël Wiertlewski

Interferometric tribometer for wide-range/high-bandwidth measurement of tactile force interaction

Force sensors are important instruments for studying the interaction between users and stimuli, and gain insight in tactile perception. This paper presents a high-precision, high-bandwidth force sensor that relies on an interferometer to achieve millinewton precision over a frequency bandwidth that spans continuous forces all the way up to 1300 Hz. The picometer resolution of the interferometer allows to design a stiff elastic structure while preserving a sufficient force measurement sensitivity.

I. INTRODUCTION

Recording the interaction between a person's finger and a tactile stimulus is paramount to understand the mechanics that underlie tactile perception. Cutaneous tactile perception responds to skin deformation that is caused by the forces at the interface between the skin and the environment. These forces carry important information about the nature of the interaction. The low frequency spectrum (DC to ≈50 Hz) carries information about the perception of shape [START_REF] Janko | On frictional forces between the finger and a textured surface during active touch[END_REF], whereas the high-frequency end of the force spectrum (≈50-800 Hz) mediates textures and transient events [START_REF] Wiertlewski | On the 1/f noise and noninteger harmonic decay of the interaction of a finger sliding on flat and sinusoidal surfaces[END_REF], [START_REF] Grigorii | High-bandwidth tribometry as a means of recording natural textures[END_REF]. A common method involves recording lateral and normal forces applied on the surface to access to the friction coefficient variations. These force measurements are also essential for the design of virtual textures on surface-haptic devices [START_REF] Vardar | Effect of waveform on tactile perception by electrovibration displayed on touch screens[END_REF], [START_REF] Messaoud | Textile fabrics texture: from multi-level feature extraction to tactile simulation[END_REF], [START_REF] Huloux | Overcoming the variability of fingertip friction with surface-haptic force-feedback[END_REF].

A popular method for measuring the interaction uses strain-gauge. These sensors can capture the evolution of the frictional force on different surfaces [START_REF] Smith | Deployment of fingertip forces in tactile exploration[END_REF]. However, the elastic body on which the strain-gauge are mounted needs to be compliant enough to provide a sufficient sensitivity yet rigid enough so it does not resonate at a frequency of interest.

To study subtle fluctuations and transients, strain gauges sensors are often substituted by piezo-electric sensors [START_REF] Wiertlewski | On the 1/f noise and noninteger harmonic decay of the interaction of a finger sliding on flat and sinusoidal surfaces[END_REF], [START_REF] Shultz | The application of tactile, audible, and ultrasonic forces to human fingertips using broadband electroadhesion[END_REF]. These sensors offer a stiffness in the order of 10 6 N/m, with sensitivity in the order of 1 mN. However the inherent transduction principle makes them ill-suited to measure forces with frequencies lower than a few Hertz.

Here we present an alternative approach based on interferometry, which is unaffected by the aforementioned trade-off and is able to measure the reaction force of a finger onto a stimulus that can weight up to 200 g with a noise floor lower than a milliNewton, while at the same time providing a frequency bandwidth higher than 1 kHz all the way down to continuous forces. At its core, any force sensor is composed of a transducer that converts the deformation of the elastic body into an electrical signal. The design of a force sensor has to consider the performance of each elements.

The sensor shown Fig. 1 aims at being superior in sensitivity and frequency bandwidth to what a human is able to perceive [START_REF] Verrillo | Sensation magnitude of vibrotactile stimuli[END_REF]. It is optimized to cover a frequency bandwidth that spans continuous forces to stimulation up to f m = 800 Hz, while at the same time being able to resolve forces which amplitudes are lower than 1 mN. The desired frequency response of the system leads to a constrain on the dynamics of the system. The studied samples have a typical mass of m = 200 g which imposes the stiffness of the flexures to be in the order of k = (2π f m )2 m = 5.10 6 N/m to realize a structure that resonate at high enough frequency not to color the signal. The required stiffness is comparable to that of the piezoelectric force sensors and hundred-fold higher than 6-axis force sensors based on strain-gauges with similar sensitivity.

A stiffer structure will trade-off force sensitivity. The sensor is designed to resolve a force of ∆ f = 1 mN, and therefore the transducer has to be able to resolve a displacement of ∆x = k∆ f = 200 pm. The interferometry principle used in this sensors is one of the few transduction principle able to resolve this level of displacement.

B. Mechanical structure

The studied texture (or a glass plate) is fixed to the top of an aluminum support. The support is a lightweight yet stiff honeycomb structure. The aluminum support is suspended by three flexures that allow for lateral displacement but are virtually infinitely stiff in the other directions. This structure guides the deformation along the sensing axis of the transducer. The lateral force sensor is then suspended to the ground with two four-bar flexure linkages made out of brass on each side. Similarly to the lateral flexures, these horizontal brass flexures allow for a vertical movement while limiting motion in the other degrees of freedom.

C. Deformation measurement

The deformation of the structure is measured by a 3-axis Fabry-Perot interferometer (IDS 3010, Attocube, Munchen, Germany). The laser is guided by an optical fiber and focused by a lens (D4/F1, Attocube, Munchen, Germany), which is mounted on the base of the structure. The laser is then reflected to the lens by a mirror, mounted on the flexible part of the structure. Micrometric adjustment screws help align the laser path. The picometer-sized displacements of the three axis (2 normal, one tangential to the fingertip exploration) are then converted into analog output refreshed at 10 MHz that can be collected by a digital acquisition board.

D. Calibration

The calibration was done by first applying known force on the upper part of the sensor with a mass and then use the same weight with a string and pulley system to apply a tangential load. The calibration also numerically removes the remaining mechanical cross-talk (around 2%).

III. CHARACTERIZATION

A. Frequency Response

The frequency response of the structure to a tangential impulse is shown Fig. 2. The upper support was impacted with a hammer. The frequency response shows a main cutoff frequency around 1300 Hz, with some normal modes between 300 Hz and 600 Hz. The response is repeatable and numerical method can be apply to flatten the response over the DC-1300 frequency range.

B. Sensitivity

The measured noise floor of the presented structure and interferometer system is 0.4 mN. The range of the interferometer allows for forces above 100 N but for realistic approach the sensor will be used to measure forces with the 1 N range. In this context it corresponds to a signal to noise ratio (SNR) equal to 3.8 • 10 6 which is superior to most piezo-electric force sensors [START_REF] Wiertlewski | On the 1/f noise and noninteger harmonic decay of the interaction of a finger sliding on flat and sinusoidal surfaces[END_REF].

IV. EXAMPLE OF TEXTURE MEASUREMENT Fig. 3 presents a typical measurement of lateral and normal forces during the tactile exploration of a virtual rendering of a sinusoidal texture displayed on an ultrasonic surfacehaptic device. The subject's finger goes back and forth along a glass plate fixed on the sensor upper support. The glass plate vibrates with an ultrasonic carrier modulated by a 20 Hz sine wave. We selected this example because the modulation of the frictional force is subliminal but the sensor is able to clearly record its fluctuations.

V. DISCUSSION AND CONCLUSION

We present a custom-built 3-axis force sensor based on a Fabry-Perot interferometer able to measure continuous forces as well as fluctuations up to 1300 Hz with a better sensitivity (SNR) than the strain-gauge based solution and comparable with piezoelectric force sensor, which are insensitive to continuous forces. Better design of the elastic structure can avoid the unwanted resonant mode that affects the frequency response. This device will be a valuable instrument to analyze a wide range of natural textures and materials with the hope of being able to realistically recreate virtual counterparts on haptic surfaces.

Fig. 1 .

 1 Fig. 1. Picture of the tribometer. The sample is mounted on a honeycomb structure that provides a light and stiff support. The test sample is suspended by a set of three four-bar flexures which deflection is measured with an interferometer via fiber optics.

Fig. 2 .

 2 Fig. 2. Frequency response of the sensor. The median of ten normalized lateral displacement amplitudes is printed in black and its first and third quartiles are shaded in grey. Red dashed lines show the ±3 dB range.

Fig. 3 .

 3 Fig. 3. Haptic signal -sine texture recreated with ultrasonic friction levitation -measured with the force sensor. The lateral force is printed in blue and the total normal force in red.

Aix Marseille Univ, CNRS, PRISM surname@prism.cnrs.fr

PSA Groupe name.surname@mpsa.com

Aix Marseille Univ, CNRS, ISM

TU Delft, The Netherlands m.wiertlewski@tudelft.nl II. DESCRIPTION OF THE FORCE SENSOR A. Design rationale


Abdurrahman Erciyas

Fedai Inanir

AC Loss Calculation of Synchronous Generators Made of HTS Superconducting Armature Windings and Permanent Magnet Rotor

In this study, a novel FEM calculation method based on A-Formulation of Maxwell Equations was proposed to evaluate the transient AC current losses of superconducting generators. The calculations are carried out in 2-Dimensions. The model generator under consideration is a three-phase synchronous generator with a 12-pole and 12-slot rotor made of permanent magnets and the stator windings made of high-temperature superconducting coils. Superconducting coils are wound with 2G HTS coated conductors and the critical current intensity is around 300 A at T = 77 K. Rotor magnets generate magnetic fields in the xy-plane and only current is induced in the z-direction inside the HTS coils. The Loss calculation of the superconducting generator were then repeated with the H-formulation and the A-H Hybrid model. AC loss was calculated for different rotor speeds using all three calculation approaches. The induced current distributions and magnetic field distributions in HTS coils for various time instants for 60 rpm rotor sped are presented.

Introduction

There is a general trend that superconducting rotationary machines of small sizes and high power are more viable commercially. The main reason is that the cryogenic subsystem required to be added to the system constitutes a much smaller part of the total machine weight and cost. The commercial availability of stranded NbTi superconductors since the late 1960s made it possible to start applications for superconducting rotating machines (Woodson et al., 1966). Since these superconductors are only suitable for carrying direct current (DC), they were originally designed to replace the DC field windings in existing rotating machines. The ideal applications for these designs were large synchronous motors, generators and unipolar DC machines.

The magnetic field strength in the air gap of conventional rotary machines is generally limited by the saturation of the iron core teeth on the rotor and stator. However, the high current carrying capacity of superconductors made it possible to generate much higher air gap field strength with a small volume field winding and very little energy loss. First, the resistance losses are equal to the field and armature windings. From this, the promised efficiency improvement can be achieved by replacing the superconducting windings with the copper field windings (representing about half of the total copper losses). In this sense, large synchronous generator designs are General Elekcric (Keim, 1985), Westinghouse (1977), Super-GM (Nitta, 2002;Rogalla & Kes, 2011) etc. has been implemented by companies. Although Super-GM implemented the project in the 1990s in which it tested three different rotors 70 MVA generators using three different NbTi field windings, none of them were found economically attractive due to technical problems such as inadequacy and cost of low temperature cooling systems, instability of the superconductor wire.

The discovery of HTS in 1986 provided an impetus for the development of various electrical machines based on these new superconductors. In particular, the fact that these superconductors offer the opportunity to operate at higher temperatures (between 25 K and 77 K) with lower cost cooling systems has enabled applications where YSS technologies such as full-scale high or low speed motors and generators are used worldwide.

Framework of Modelling

New calculation methods have been developed for electromagnetic analysis of superconductor generators. The calculations were made in two dimensions and the model generator was designed. Figure 1. The model design considered is a three-phase synchronous generator with 12 poles and 12 slots made of permanent magnets and the stator windings made of superconducting coils. Superconducting coils are 2G coated conductors wound and the critical current intensity is around 300 A at T = 77 K. Rotor magnets generate magnetic fields in the xy-plane and only current is induced in the z-direction from the HTS coils. Parameters used for model generators are given in Table 1. Firstly, A-V formulation was used for calculations. According to the A-V formulation, the current induced in the generator is calculated according to the Ampere equation as follows:

𝜎 𝜕𝑨 𝜕𝑡 + 1 𝜇 0 𝜇 𝑟 ∇ × ∇ × 𝑨 = 0 (1)
where σ represents the conductivity of the superconductor and the materials used with the generator, A represents the vector potential, μ 0 the magnetic permeability of the cavity, and μ r the relative magnetic permeability. A new method was developed to solve the time dependent partial differential equation with the finite element software COMSOL Multiphysics software.

Figure 1. Model generator used in calculations. The rotor of the generator is designed from permanent magnets and the magnets are positioned to generate magnetic fields in the radial direction inward and outward.

Fedai Inanir acknowledges the financial support of the Scientific and Technological Research Council of Turkey (TUBITAK) under the contract number 114F424.


Shuangrong You
email: shuangrong.you@vuw.ac.nz

Swarn S Kalsi
email: swarnkalsi@yahoo.com

Mark D Ainslie
email: mark.ainslie@eng.cam.ac.uk

Rodney A Badcock

Nicolas J Long

Zhenan Jiang

AC Loss Simulation in HTS Armature Windings of A 100 kW All-HTS Motor

Keywords: HTS motor, AC loss, T-A formulation, moving mesh

ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

High AC loss in armature windings is one of key issues for achieving all-HTS motors [START_REF] Jiang | Numerical analysis and finite element modelling of an HTS synchronous motor[END_REF]. As a powerful tool, FEM modelling has been widely used to study the electromagnetic behaviours in HTS motors in recent years [START_REF] Vargas-Llanos C R | T-A formulation for the design and AC loss calculation of a superconducting generator for a 10 MW wind turbine[END_REF][START_REF] Zhang | AC loss estimation of HTS armature windings for electric machines[END_REF][START_REF] Zhang | Design and Modeling of 2G HTS Armature Winding for Electric Aircraft Propulsion Applications[END_REF]. In this work, a 100 kW, 1500 rpm all-HTS pump motor was numerically studied in order to investigate the loss behaviours in the armature windings.

The schematic of cross-section in this motor is shown in the figure 1. The main design parameters for this HTS motor are listed in Table 1.

Numerical modelling and results

In this motor, both field windings and armature windings are designed using REBCO conductors. However, the AC loss in the field windings is not considered. The model was simplified by only applying E-J power law and asymmetric Ic() characteristics of the HTS wires in the armature windings while the current density in the field windings was set as a constant. The T-A formulation was applied to this motor model [6 -8]. T formulation was applied to the HTS armature windings by using "Coefficient Form Boundary PDE" interface in COMSOL and A formulation was applied to whole FEM model by using "Rotating Machinery, Magnetic" interface. In order to study the electromagnetic behaviours of whole motor at a rotating condition, moving mesh was applied to the rotor part of the model.

In our current work, 14-turn double race-track coils wound with 4 mm wide REBCO wires were considered as the armature windings. The distribution of current density normalized by the local 𝐽 𝑐 in one slot at 1/4 circle is shown in figure 2. The result with current is slightly greater than that without current. This implies most of AC loss in the armature windings is magnetization loss. Transport AC loss in the windings is negligible because the designed operating currents in an HTS motor is much lower than its critical current ( 𝐼 𝑡 𝐼 𝑐 is smaller than 0.5) for safety consideration. Various coil designs will be considered in the following works to reduce AC loss in the HTS armature coil windings.

Figure 1 .

 1 Figure 1. Schematic cross section of the 100 kW, 1500 rpm all-HTS motor

  The asymmetric Ic() characteristics of the HTS wires were from the Robinson Research Institute HTS wire data base [5].

Figure 2 .

 2 Figure 2. Distribution of normalized current density in one slot of armature windings at 1/4 circle. Each coil has 7 × 2 = 14 turns HTS wires.Magnetic field distributions of the entire motor are illustrated in figure3. Due to the iron poles yoke, HTS windings can avoid experiencing high magnetic flux density, which in turn reduces the AC loss in the HTS windings.

Figure 3 .

 3 Figure 3. (a) at 1/6 and 1/3 circle, distribution of magnetic flux density norm in whole motor.The instantaneous AC loss for all armature windings with/without operating current is shown in figure4.

Figure 4 .

 4 Figure 4. (a) Instantaneous loss in armature windings with/without operating current (1.5 kA ampere-turn).

Table 1 .

 1 Main parameters of the HTS motor Ic of the HTS wire at 65 K s.f. (A/cm) 765

	Parameters	Value
	Output power (kW)	100
	Rotating speed (rpm)	1500
	Operating temperature (K)	65
	Ampere-turns per armature coil (kA)	1.5
	Ampere turns per field coil (kA)	4.0
	Outer diameter of the motor (mm)	200
	Air gap (mm)	2
	Diameter of the rotor (mm)	120
		Silicon Steel
	Material of iron core and yoke	Non-grain
		oriented M-36


Keyang Wang

V A Anvar

Arend Nijhuis

Yuanwen Gao

Jeremy Weiss

Danko C Van Der Laan

St K Wang

M S A Hossain

Analysis of contact behavior in CORC® cabling and under axial tension

Keywords: CORC ®, cable, wire, REBCO, HTS, Finite Element Method, modeling, normal contact force, friction, critical current

à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

INTRODUCTION

The work on CORC ® cable and wires is categorized into two sections.

A. Winding process of CORC ® and Pancake coil

To produce CORC ® cables and Pancake coils, the REBCO tapes are wound around a center core with a particular helix angle or under only few degrees. During the winding process, the tape's deformation is affected remarkably by the contact with the core. This affects the critical performance because the REBCO tape is strain sensitive. The mechanical properties of Pancake coils and CORC ® cables are related to the initial contact during winding. The interlayer contact resistance is also a critical electrical property, since it determines the ability of current sharing and therefor the cable's stability. However, the winding is a rather complex process because of the combined deformations and its strong non-linearity, which leads to difficulties in calculating the contact force.

A two-dimensional simplified winding model of a Pancake coil is established through theoretical methods, considering the two main factors in the winding process; Relaxation and Poisson effect (longitudinal stress produces transverse deformation). The model contains two defined main regions: the core and the REBCO tape. The core is the rigid domain and the tape is the solution domain. One end of the tape is fixed to the core's surface and the other end is given a pre-tension force Fpre. The rigid core is designated to rotate. The model is shown in Fig 1. As a result, the contact behaviour between the tape and the core is obtained. Then, the three-dimensional winding model is established by COMSOL software to analyse the contact behaviour between HTS tapes and center core. After describing the Pancake coiling and the CORC ® cabling processes, the contact stress distribution is calculated using the FEM model (see The tape used in the test is SCS4050 (4 mm width and 50 m substrate thickness), which is provided by SuperPower Inc. The experimental test is performed by hanging different weights (50,100,200,300,500,700,900,1000,1500,2000,2500, 3000 and 4000 g) acting as the winding pre-tension force, and then slowly rotating the aluminium roller (r=10 mm). The measured result shows that as the winding pre-tension force increases, the Relaxation effect becomes weaker. The results obtained from the theoretical calculations, numerical simulations and winding tests are consistent. However, the numerical model is smaller than the theoretical calculation result and closer to the test data because of the plasticity. B. Tension process of CORC ® cables.

After the cabling process, the CORC ® cables or wires will be in further processing stages, such as stretching, bending, and winding. The HTS tapes in the cables will be deformed by interlayer interaction, resulting in possible performance degradation. The contact behaviour, both the tape-to-tape contact and the tape-to-core contact in cables and coils, are more complex. The study of this contact behaviour not only helps in understanding the mechanics of the cabling process but also provides the basis for studying the current distribution in cables and coils. The FEM results show the normal contact stress and the distribution of the REBCO strain for each layer under a total tensile strain of 0.5, 1.0, 1.5, and 2.0%. (The average normal contact force of each layer is shown in Fig 6). Firstly, the normal contact stress of each layer increases as the tensile strain increases. The inner layer is subjected to more normal contact pressure than the outer layer, and the innermost layer is the largest. This is because when stretching the CORC ® wire, the innermost layer withstands the cumulative sum of all the layers above it.

Fig 1 :

 1 Fig 1: Pancake coil and CORC ® cable winding model with boundary conditions.

  Fig 2).

Fig 2 :

 2 Fig 2: The process of CORC ® cabling.(r=3 mm, T=100 MPa, 2w=4 mm, t=0.096 mm, α=45 o ). Then the FEM results are compared with the theoretical results and verified by experimental winding tests. A test platform, as shown in Fig 3, is designed and produced to verify the theoretical and numerical simulations. Fig 4 compares the relaxation offset angle and relaxation offset distance under different winding pre-tension force through theoretical calculations, numerical simulations and experimental tests.

Fig 3 :

 3 Fig 3: (a) The REBCO tapes produced by SuperPower Inc. (b) The used weights. (c) Relaxation effect test. (d) Schematic diagram of measurement method.

Fig 4 :

 4 Fig 4: (a) Relaxation angle θ. (b) Relaxation offset distance b with different winding pre-tension stress.

Fig 5 .

 5 Fig 5. The boundary condition of the FEM model. In order to analyze the contact behaviour of CORC ® wires under axial tension, a three-dimensional numerical model having 30 tapes arranged in 12 layers is developed, as shown in Fig 5. A simple theoretical model can estimate the tape strain, the normal contact and friction forces during tensile loading, which helps to understand the layer interactions.


Ze Jing
email: zejing@xidian.edu.cn

Mark Ainslie
email: mark.ainslie@eng.cam.ac.uk

Coupled multiphysics modelling of the thermalmagnetic-mechanical instability in bulk superconductors during magnetization

Keywords: multiphysical modelling, bulk superconductor, magnetization, mechanical failure, pulsed field magnetization I

Magnetization is one of the most significant challenges for the application of bulk superconductors. During the pulsed field magnetization process, large numbers of quantized magnetic flux vortices can rush in or out of the superconducting sample, which generates a significant amount of heat and can lead to the thermomagnetic instability of the bulk superconductor. The large Lorentz force and the sharp temperature rise can result in mechanical failure (or instability), which also hinders the application of bulk superconductor. In addition, thermomagnetic instability and mechanical failure may couple with each other. Thus, it is challenging to simulate the magnetization of bulk superconductors used as high field magnets. In this paper, a numerical simulation framework based on the coupled H-formulation for the electromagnetic behavior of superconductors, the heat diffusion equation, and the phase-field model for the failure of solids is proposed and implemented to simulate the thermal-magnetic-mechanical instability behavior of bulk superconductors during the magnetization process. The thermal-magnetic-mechanical instability-induced dendritic flux patterns are obtained and the corresponding material damage are evaluated through these numerical simulations. The magnetic field, current density, temperature, stress/strain field distribution, and the material failure within the bulk sample in the high magnetic field magnetization process will also be presented.

INTRODUCTION

Large, single-grain REBa 2 Cu 3 O 7-δ (REBCO) bulk high temperature superconductors, which can trap magnetic fields up to 17.6 T at 26 K [START_REF] Durrell | A trapped field of 17.6 T in melt-processed, bulk Gd-Ba-Cu-O reinforced with shrink-fit steel[END_REF], have demonstrated significant potential for use in a variety of engineering applications [START_REF] Hull | Applications of bulk high-temperature Superconductors Proceed[END_REF]. However, REBCO is a ceramic material, brittle in nature with poor mechanical properties [START_REF] Krabbes | High temperature superconductor bulk materials: fundamentals, processing[END_REF]. Mechanical stress is inevitable for their practical use, which needs special attention especially for high magnetic field applications. The range of applications of bulk superconductors is limited by their poor mechanical stability. Fracture and microscopic damage have been found extensively in REBCO bulk superconductor magnets during the high magnetic field magnetization process. Moreover, the motion of magnetic flux during the magnetization process dissipates a lot of energy and can lead to thermomagnetic instability. Thermomagnetic instability induced flux jumps are frequently observed during pulsed field magnetization.

Giant magnetostriction (up to 10 -4 ) has been measured by Ikuta et al [START_REF] Ikuta | Giant magnetostriction in Bi2Sr2CaCu2O8 single crystal in the superconducting state and its mechanism[END_REF]. Cracking of a mini-magnet composed of four YBCO discs under an applied field of 14 T at 49 K has been observed by Ren et al [START_REF] Ren | Damage caused by magnetic pressure at high trapped field in quasi-permanent magnets composed of melt-textured YBaCuO superconductor[END_REF]. Takahashi et al [START_REF] Takahashi | Fracture behavior analysis of EuBaCuO superconducting ring bulk reinforced by a stainless steel ring during field-cooled magnetization[END_REF] observed fracture behavior of an EuBaCuO ring bulk during field-cooled magnetization and numerically simulated the mechanical stress. Ainslie et al developed a modelling framework that couples electromagnetic, thermal and structural mechanics, and analyzed the mechanical stresses and flux jumps of bulk superconductor during magnetization [START_REF] Huang | Numerical modelling of mechanical stresses in bulk superconductor magnets with and without mechanical reinforcement[END_REF][START_REF] Fujishiro | Trapped field of 1.1 T without flux jumps in an MgB2 bulk during pulsed field magnetization using a split coil with a soft iron yoke[END_REF]. Research group lead by Zhou has made significant contributions to modelling the mechanical behavior of bulk superconductors [START_REF] Wu | Stress analysis in high-temperature superconductors under pulsed field magnetization[END_REF], and has numerically simulated the fracture [START_REF] Chen | XFEM analysis of the fracture behavior of bulk superconductor in high magnetic field[END_REF] and flux jump [START_REF] Xia | Numerical investigations on the characteristics of thermomagnetic instability in MgB2 bulks[END_REF] behavior of the bulk superconductor.

In this paper, a theoretical modelling framework which couples the H-formulation, the heat diffusion equation, and the phase-field fracture model is proposed to simulate the coupled electromagnetic-thermal-mechanical behavior of bulk superconductors during the magnetization process.

II. THEORETICAL MODELING AND NUMERICAL METHODS

A. Electromagnetic modelling

The electromagnetic behavior of the bulk superconductor is governed by the following Maxwell's equations:

∇× = H J
(1)

t ∇× = -∂ ∂ E B /
(2) where H is the magnetic field strength, J represents the current density, E is the electric field. The magnetic flux density B is related to H by B=μ 0 H, where μ 0 is the vacuum permeability.

The relationship between the electric field and the current density is given by the E-J power law,

1 0 n c c E J J - = E J / J / (3)
where E 0 is the electric field constant (E 0 =1×10 -4 V/m), J c is the critical current density and defined as a function of the temperature T, the magnetic flux density B and the phase field variable φ by:

( ) ( ) ( ) ( ) 2 1 2 m a x m a x 1 ( , , ) 1 exp / / exp 1 / y c c Lc J B T J B B J B B B B y φ φ     = - - + -         (4)
in which J c1 , J c2 , B L , B max and y are temperature dependent parameters and can be found in [START_REF] Huang | Numerical modelling of mechanical stresses in bulk superconductor magnets with and without mechanical reinforcement[END_REF] and references therein.

B. Heat diffusion equation

The thermal behavior is described by

( ) ( ) 0 / c T t T h T T κ ∂ ∂ =∇⋅ ∇ - -+ ⋅ J E ( 5
)
where κ is the thermal conductivity, c is the specific heat, T 0 is the ambient temperature and h is the heat transfer coefficient.

C. Phase-field modelling

In the phase-field fracture model, the coupled governing equations for displacement and the phase field are:

( ) 0 0 / 2 / 2 c c G l Gl in φ φ + - Δ = Ω   (7)
The boundary conditions along the domain are given by: t on ⋅ = ∂Ω σ n t [START_REF] Fujishiro | Trapped field of 1.1 T without flux jumps in an MgB2 bulk during pulsed field magnetization using a split coil with a soft iron yoke[END_REF] u on = ∂Ω u u [START_REF] Wu | Stress analysis in high-temperature superconductors under pulsed field magnetization[END_REF] 0 on φ ∇ ⋅ = Γ n [START_REF] Chen | XFEM analysis of the fracture behavior of bulk superconductor in high magnetic field[END_REF] where n denotes the outward normal vector to the boundary, and u and t are the prescribed displacement and traction force on the boundary, respectively.

In our numerical simulations, a multiphysical field coupling scheme based on the finite element method is proposed and applied in the commercial software COMSOL Multiphysics. The electromagnetic equations are solved using a general form PDE, the phase field fracture model is adapted into the Solid Mechanics module and the heat equation is applied in the Heat Transfer module.

III. RESULTS AND DISCUSSIONS

Here, we present the numerical simulation results for the thermal-magnetic-mechanical instability behavior of a cylindrical bulk superconductors with reinforcement during pulsed field magnetization. The variations of electromagnetic fields in the thickness direction is neglected, and a 2D planestrain problem is considered in our simulation. The radius of the bulk sample is 10 mm, and thickness of the stainless-steel reinforcement ring is 2 mm. The practical experimental waveform for the applied external magnetic field B app observed during the PFM is taken as:

 exp 1  .
We take τ = 12 ms in the simulation, typical of experiments. Fig. 1 shows the time evolution of magnetic flux (a 1 )-(d 1 ), the temperature (a 2 )-(d 2 ), the phase field (a 3 )-(d 3 ), and the critical current density (a 4 )-(d 4 ) during a typical pulsed field magnetization process with B max = 7 T at 40 K. It is clearly shown that dendritic flux avalanches (similar to that found in thin disks [START_REF] Vestgården | Lightning in superconductors[END_REF]) develop within the bulk sample from the periphery as the applied field rises to about 4 T at 40 K. In Fig. 2, we present the magnetization jump curve and the maximum temperature of the bulk superconductor during pulsed field magnetization. It is shown that the magnetization changes sharply during the avalanche process, and the temperature can even rise up to 180 K, which may induce mechanical damage in the bulk (as shown in Fig. 1

Fig. 1 .

 1 Fig. 1. Magnetic flux, temperature, phase field and the critical current density distribution within the bulk superconductor with reinforcement under pulsed field magnetization (B max = 7 T at 40 K).

  (a 3 )-(d 3 )).

Fig. 2 .

 2 Fig. 2. Magnetization jump and the maximum temperature of the bulk superconductor under pulsed field magnetization with maximum magnetic field B max = 7 T at 40 K.

ACKNOWLEDGMENT ZJ acknowledge the support from the National Natural Science Foundation of China (Nos. 11602185 and 11972271), Natural Science Basic Research Plan in Shaanxi Province of China (Grant No. 2018JQ1013) and the financial support from China Scholarship Council. MA would like to acknowledge financial support from an Engineering and Physical Sciences Research Council (EPSRC) Early Career Fellowship, EP/P020313/1. All data are provided in full in the results section of this paper.


Abderrahmane Djaouti
email: abderrahmane.djaouti@univ-lorraine.fr

Gaëtan Didier
email: gaetan.didier@univ-lorraine.fr

Denis Netter
email: denis.netter@univ-lorraine.fr

Modelling and simulation of resistive type superconducting fault current limiter using OpenModelica

Keywords: Superconducting fault current limiter, hot spot, finite difference, OpenModelica, electro-thermal model

Modelling and simulation of resistive type superconducting fault current limiter using

OpenModelica.

I. INTRODUCTION

The growing demand for electrical energy leads to increasing of current levels, which can exceed the maximum permissible values of the electrical equipment. As a consequence, the high thermal and mechanical stresses caused by fault currents can damage electrical devices. The new generation of SFCL, called the second generation 2G HTS is one of the most common ways to limit high-fault currents.

In fact, the FCL is essentially a variable impedance that is connected In series with a circuit breaker. During a fault, the impedance value increases until the fault current is reduced to a lower level that the circuit breaker will be able to handle. In addition to short-circuit protection, the application of FCL in electrical power systems can improve the stability of the system in the case of disturbances .

Since the discovery of high-temperature superconductors (HTS) in 1986, the superconducting fault current limiter has been considered as the most promising and innovative application in power systems. The main aspect to develop SFCL systems is the understanding the behavior of superconducting materials in quenching transition. At present, different rSFCL models have been developed. However, the majority of them are limited to a purely electrical model or a 1D electrothermal model that takes into account a homogeneous temperature distribution over the thickness or length. This paper proposes a modeling and simulation of high temperature resistive SFCLs using OpenModelica software considering the electrical and thermal behavior along the length and thickness (2D). The heat equation is discretized by the finite difference method, considering the heat propagation according to the length and thickness of the HTS. This equation combined with the eclectic algebraic equations represents a system of differential/algebraic equations (DAE). The solving method : differential/algebraic system solver (DASSL) integrated in OpenModelica represents a main advantage to use a such software. OpenModelica offers also the advantage that it is compatible with standard Functional Mock-up Interface (FMI) for co-simulation. FMI is a toolindependent standard for model exchange and co-simulation.

II. MODELING 2D OF HIGH TEMPERATURE SUPERCONDUCTORS TAPE

A. Description of the HTS tape

The model of HTS tape contains a stack of layers: buffers (oxides) and a superconducting layer e.g. (RE)BCO which is deposited on a metal substrate (stainless steel and Hastelloy). The entire tape is coated with a silver layer. However, in order to ensure a two-dimensional (2D) aspect, no silver layer has been added to the lateral side of the tape. Fig. 1 presents the discretization used for the HTS tape considering the interfacial resistance (𝑅 𝑖𝑛𝑡𝑒𝑟 ) between the superconducting and silver layers. The length and thickness of each element depends on the number of discretization selected in the thickness and in the length.

Silver

B. Electrical Model

The resistance of each discretization element of HTS tape can be modeled as a non-linear resistance (1) depending on the current I and/or the temperature T.

𝑅 𝑒𝑙 (𝐽, 𝑇) = 𝜌(𝐽, 𝑇) 𝐿 𝑆 (1)

Where 𝜌(𝐽, 𝑇) , L and S are respectively, the resistivity, the length and the longitudinal cross-section of the considered element. The resistivity of the superconducting layer has been evaluated by the non-linear E-J characteristic called the power law [START_REF] Bonnard | Multi-scale model of resistive superconducting fault current limiters based on 2G HTS coated conductor[END_REF]. The electrical algebraic equations are obtained by applying Kirchhoff's circuit laws.

C. Thermal Model

The diffusion of heat both in thickness and length is important, especially when the fault current is close to the critical current 𝐼𝑐 or in the case of hot spot phenomena [START_REF] Bonnard | Multi-scale model of resistive superconducting fault current limiters based on 2G HTS coated conductor[END_REF]. The model proposed in this paper takes into account the twodimensional (2D) behavior of heat conduction in the thickness and length (2D) of the HTS strip. The partial derivatives of the heat equation (2) have been transformed to a finite difference discretized form.

where 𝑄𝑗 represents the total Joule losses inside the considered element, 𝜅 is the heat conductivity of the element, 𝜌𝑚 is the mass density and 𝐶𝑝 is the specific heat capacity.

III. SIMULATION RESULTS

The presence of critical current density 𝐽 𝑐 inhomogeneities in superconductors causes the most critical defect (hot spot). To validate the proposed model, the power grid and the parameters studied in [START_REF] Bonnard | Multi-scale model of resistive superconducting fault current limiters based on 2G HTS coated conductor[END_REF] have been used. A tape of 5 mm in length having a hot spot (𝐽 𝑐𝑓 = 0.55𝐽 𝑐 ) of 1 mm in length (Fig. 2) supplied by a DC source equal to 100A is used.

Silver

Fig. 1 .

 1 Fig. 1. HTS tape discretization.

  𝜌 𝑚 𝐶 𝑝 (𝑇) 𝜕𝑇 𝜕𝑡 = 𝛻. (𝑘(𝑇)𝛻𝑇) + 𝑄 𝑗 (𝑇, 𝐼)

TABLE I .

 I PARAMETERS OF THE (RE)BCO SUPERCONDUCTOR[START_REF] Bonnard | Multi-scale model of resistive superconducting fault current limiters based on 2G HTS coated conductor[END_REF]

	Parameters	Value	Description
	Ec	1 V cm -1	Critical electrical field criterion
	Jc0	2.5 MA cm -2	Self-field critical current density
	n	15	Power law exponent
	𝜌𝑇c	30  cm	Normal state resistance at 𝑇c
	α	0.47  cm -1	Temperature coefficient
	𝑇c	90 K	Critical temperature
	T0	77 K	Temperature of the LN2 bath
	TMax	250 K	Maximum admissible temperature

The solver (DASSL) integrated in OpenModelica has been used to solve the non-linear coupled equations with a large number of variables.