Nicolas Lamaison

David Chèze

Florent Lefrançois

Fabien Bruyat

Jean-François Robin

Operational Behavior of a Solar-Fed Bidirectional Substation for 4GDH Networks

Keywords: Solar Feed-In, Decentralized, Bidirectional Substation, District Heating

Decentralized feed-in of solar surplus heat from prosumers into a District Heating Network is here addressed. A specified architecture combining a Return/Supply connection, a set of valves, a feed-in pump and a primary side connection of the solar field led to the conception and fabrication of a bidirectional substation prototype. The latter is tested in a hardware-in-the-loop configuration regarding the consumer and solar field sides while the DHN side of the substation prototype is connected to a real DHN. A 12-days sequence, representative of typical days, with varying solar irradiation, DHN thermo-hydraulic conditions and consumer demands is performed with the prototype. A set of indicators showed i) an inappropriate behavior of the controller during the coldest days, ii) various improvements of the solar fraction during the hottest days (with up to 52% of total solar fraction) and iii) rather good performances of the controller (about 90% effectiveness). The analysis was completed by the observation of detailed dynamic result for one day. Different operating regimes (full DHN consumption, solar self-consumption complemented by the DHN, solar self-consumption and feed-in) could be observed during a single day highlighting the modularity of the BSST to switch automatically (by hydraulic balancing) from one mode to the other.

Introduction and State-Of-the-Art

In the "2way District Heating" course of action from the 4th Generation District Heating concept [START_REF] Lund | 4th Generation District Heating (4GDH)[END_REF], decentralized feed-in of solar heat from prosumers seems to be a promising solution to increase the share of renewable energy in District Heating Networks (DHN), especially in dense urban areas with limited ground surface and high cost of land. However, when scattered customers roofs are used to collect and inject heat locally in a network, a new component is required. Such a component, referred here as Bidirectional SubSTation (BSST), must be able to transfer heat i) from the DHN to the local consumer when the locally produced solar heat is not sufficient and ii) from the local consumer to the DHN when solar heat is in excess.

Among the various reinjection principles, Return to Supply (R/S) feed-in is chosen since it seems to be the most flexible option from the DHN point of view [START_REF] Lennermo | Distributed heat generation in a district heating system[END_REF]. Indeed, it allows a reasonably low temperature in the solar field (better panels' efficiency) without increasing the return line temperature (better efficiency of other heat generators and lower heat loss) while presenting lower risks in terms of pipe fatigue. However, R/S feed-in implies the necessity to cope with the variable i) local differential pressure between the return and supply lines, ii) solar irradiation, and iii) local DHN return temperature. Moreover, the feed-in temperature must be close to the local supply line temperature of the network to prevent both pipe fatigue [START_REF] Lennermo | Control of decentralised solar district heating[END_REF] and unmet customer demand.

While architecture and operational behavior of feed-only substations are addressed in the literature rather extensively [START_REF] Lennermo | Control of decentralised solar district heating[END_REF][START_REF] Rosemann | DH Networks -Concept, Construction and Measurement Results of Decentralized Feed-In Substations[END_REF][START_REF] Schäfer | Technical challenges for solar thermal plants with decentralized feed-in into district heating networks and deduced plant concept for the experimental feed-in station in the SWD[END_REF], 2-way substations studies are seldom in the open literature. Though, the combined substation concept of [START_REF] Rosemann | A New DH Control Algorithm for a Combined Supply and Feed-In Substation and Testing Through Hardware-In-The-Loop[END_REF] can be noted. Starting from usual feed-only architecture combining a set-of-valves and a feed-in pump on the primary side, they used a rather complex control algorithm and architecture of the substation combining on the secondary side the solar field together with consumer internal heating loops.

In the frame of the Horizon 2020 THERMOSS project, specifications, modeling and prototype testing of a 2-way substation with the solar field connected on the primary side is performed. While the specifications and modeling approach were presented respectively in [START_REF] Lamaison | A multi-criteria analysis of bidirectional solar district heating substation architecture[END_REF] and [START_REF] Lamaison | Modelling and Simulations of solar two-way substation[END_REF], the present paper deals with the design and original test approach of a BSST prototype combining hardware-in-the-loop and real experimental DHN.

Architecture and Control

Architecture

To build the BSST architecture, an initial review on two-way substation both from the industrial and research point of view led to the fact that i) nothing exists in the catalogue of usual manufacturers; ii) the Return/Supply (R/S) type of feed-in is the most promising compared to other types (R/R or S/S).

Discussions with DHN operator and substation manufacturer led to specifications in the form of constraints and requirements, insisting specifically on modularity and cost. Different features (local usage, network connection and control strategy) were then combined to propose various architectures. After pros and cons listing, solutions involving complete reinjection of the heat into the DHN seemed the most promising [START_REF] Lamaison | A multi-criteria analysis of bidirectional solar district heating substation architecture[END_REF], leading to the technical drawing shown in Fig. 1. Especially, the substation does not include a storage as usually associated to a solar field. Indeed, thanks to the presence of the DHN as heat sink and because of the extra heat loss induced by a storage, it is not necessary and thus preferred not to include a storage in the architecture In Fig. 1, connections to the supply and return lines of the DHN are schematically represented on the left block. Regarding the right block, it represents the consumer and exhibits a separation of space heating (SH) and domestic hot water (DHW) supplies at the primary side (dedicated heat exchangers) Regarding the middle block, it is comprised of the connection to the solar loop and the entire feed-in loop components. For the latter, the main components are the feed-in pump (P_SST_FEED), the heat exchanger (ECH_BESST_SOL) and the feed-in valve (V_BSST_FEED). Red and blue lines respectively indicate hot and cold lines. Red, blue, orange and black labels respectively denote temperature probes, pressure probes, mass flow meters and equipment. The main equipment (pumps, heat exchangers and valves), small hydraulic components (filters, check valves) and the instrumentations are highlighted. Also, it is important to mention that drains (all loops), safety valves (all loops) and expansion tanks (solar and SH loops) are not shown.

Main principles and control

The two main operating modes of the substation are schematically represented in Fig. 2:  In mode 1, the heat consumption for domestic hot water and space heating is larger than the solar production. The solar energy if any is used for the consumer needs in addition to the heat coming from the network. The flow in the service lines is thus from the supply to the return line;

 In mode 2, the heat consumption is lower than the solar production. The solar energy is used to entirely satisfy the consumer needs if any and the surplus heat is reinjected to the network. The flow in the service lines is thus from the return to the supply line.

Sizing, test setup and procedure

Sizing

Following a stage of modeling [START_REF] Lamaison | Modelling and Simulations of solar two-way substation[END_REF], the substation has been sized (heat exchangers, valves and pumps) for 42 kW in SH, 60kW in DHW and 70kW in solar. According to the sizing and architecture, the substation has been assembled as shown in Fig. The entire setup for testing the substation is finally as shown in Fig. 5 with the BSST prototype surrounded by different components.

 The emulated building space heating (SH): The building model is a typical 4 levels building erected in the early 2000's in France, parallelepiped U-shape, convex side facing south. It is composed of 12 apartments, each heated according to its thermal zone temperature by local radiator flow rate adjustment with thermostatic valve. The space heating flow is a constant value (local bypass at radiator level when thermostatic valve adjust the flow rate in radiator).  The emulated building domestic hot water (DHW) and sanitary loop: The DHW model consists in a repeated daily pattern draw-offs profile for 12 households (obtained using [START_REF] Jordan | DHWcalc: Program to generate domestic hot water profiles with statistical means for user defined conditions[END_REF], with total daily consumption variation from month to month. The DHW sanitary loop demand is not calculated from simulation but rather physically emulated. To do so, a specific hydraulic module for real heat dissipation to ambient temperature was built. The water flow rate in the loop is fixed to 300kg/h by recycling pump and the heat dissipation is controlled dynamically to keep the inlet-outlet temperature drop close to 5K.  The emulated solar field: The modeled solar field is facing south and is composed of 50 m² of common technology flat solar panel. It is worth mentioning here that the size of this solar field is limited by the capacity of the hardward-in-the-loop setup. In theory and as discussed in [START_REF] Lamaison | Modelling and Simulations of solar two-way substation[END_REF], the size of the solar field can be larger when associated to such substation.



The real DHN: Unlike the building and solar field sides that are emulated, the DHN side is connected to a small real DHN with varying local differential pressure and return temperature conditions. It is composed of a district heating network of about 200m length, supplying heat to real buildings and to an emulated building (the one in which the BSST is connected). The main production units of this DHN are a condensing gas boiler of 280kW, a heat pump of 50kW and a solar field of 300m². The network is also equipped with a hot storage tank of 40m 3 .

Results and discussion

The results are firstly presented globally for each day of the test sequence using indicators and secondly with the analysis of the dynamic behavior for a specific day of operation. When reinjection occurs, energy-based integral when setpoint is satisfied within a 2K margin divided by the total energy-based integral With the indicators described beforehand, Tab. 3 now presents those for the 12-days sequence of our testing procedure, as described in Section 3. The main conclusions to draw from this Table are the following:

 During the 6 colder days, the solar field does not produce. This is explained by the too high outlet solar field set point constrained by the reinjection arrangement. Indeed, the controller setpoint calculation is based on the assumption that reinjection may occur and is thus rather high because the network operates at high temperature (an average of 75°C for these 6 days). The latter combined with the low irradiation and the low temperature in the solar loop piping (because of the low ambient temperature) leads to an inappropriate solar field operation.

 During the 6 days for which solar reinjection takes place, we observe various levels of solar fraction but more striking various level of total solar fraction. The latter is due to the increased gap for hotter days between the local demand of the building and the energy produced by the solar field. It is worth noting here that the largest value of total solar fraction of 52.0% is still not that high. The reason is the rather small size of the emulated solar field because of test bench limitations (see Section 3). Additional tests not shown here were performed with real solar field of 70kW and solar fractions up to 80% were reached.

 Regarding the performances of the controller, the standard deviation of the feed-in temperature is rather high compared to what was expected from previous tests with a real solar field (see Section 3), for which values around 1.5°C were obtained. These results are explained by the fact that the controller parameters were calibrated with this real solar field and not refreshed after the switch to an emulated solar field (with lower aperture area and less inertia). However, the feed-in effectiveness exhibits reasonably good values (average of 92%) highlighting the fact that reinjection on the network occurs at an appropriate level of temperature.

1 0 0 98.2 0 0 0 0 0 - - 2 0 0 98.3 0 0 0 0 0 - - 3 0 0 98.4 0 0 0 0 0 - - 4

Focus on dynamic results for a specific day

In the present section, we observe the detailed dynamic results for day n°6 of the sequence with the indicators results being highlighted in Tab. 3.

Fig. 7 presents thus results obtained for day n°6 of the sequence between 9a.m. and 7p.m. in terms of various heat transfer rates and temperatures. Fig. 8 presents additional results for the same observation period in terms of heat power at DHN level, controller actuators signals and differential pressures.

Focusing on regions A, B and C highlighted in the Figure, we observe the following operating regimes:

 Region A: The BSST is in consumption mode from the DHN point of view (mode 1 in Fig. 2). The outlet solar loop temperature and thus the outlet feed-in loop temperatures are below the feed-in set point. The feedin pump is off and the solar pump operates at its lowest speed. The building consumption is thus entirely met by the DH and TT01 at the entrance of the substation on the hot service line is representative of the DHN local supply pipe temperature. The local differential pressure, here imposed by the DHN pump, is following the local building consumption indicating that i) there is no local consumption of solar heat, and ii) the other consumer on the network have a rather constant heat consumption (see top graph of Fig. 7).

 Region B: The BSST is in consumption mode from the DHN point of view (mode 1 in Fig. 2). However, periodically, the SST consumption is lower than the building consumption indicating solar heat production and self-consumption. This is the case especially at the end of the period where the outlet solar field temperature surpasses the feed-in set point. This self-consumption is achieved in a period where the substation consumption is significantly varying because of DHW draw-offs (see DHW valve opening in Fig. 8). The latter variations together with the variable solar fraction lead to large variations of the local differential pressure drop (imposed by the DHN pump) indicating strong flow rate variations in the DHN.

 Region C: The substation is in feed-in mode from the DHN point of view (mode 2 in Fig. 2). The solar heat being used for both the local building heat demand and for reinjection into the DHN. The feed-in pump is ON and the valve V_BSST_FEED together with the pump speed are adjusted to match the flow rate set point. TT09 at the entrance of the substation on the cold service line is representative of the DHN local return pipe temperature. Because of feed-in mode operation, the differential pressure is here imposed by the feed-in pump. Its variations are following the variations in demand of the other DHN consumers (see top graph on Fig. 7).

The main negative observation that can be made here for this day of operation is the low frequency periodic alternation between consumption and feed-in modes even though solar irradiation and building consumption do not exhibit such behaviors. The latter is attributed to controller setting slightly inappropriate since they were established during tests with a real solar field (with larger aperture area and more inertia) and not refreshed for this new testing setup. A very interesting insight offered by our testing setup is the observation of the DHN pressure field variations because of the reinjection. Indeed, because of the usage of a real network, furthermore linear, influence on the DHN pressure cone is directly visible. As shown in Fig. 9, which presents a qualitative explanation of the influence of the feed-in on the real DHN pressure cone, the solar heat produced at prosumer 1 benefits both the local consumer (consumer 3 = prosumer 1, cf. top graph on Fig. 7) and consumer 2, which thus reduces its dependence to the main production plant.

It is worth mentioning here that during tests with extreme conditions (larger solar field, low demand on the network), we reached situations with a fully inverted pressure cone which becomes problematic if no bypass for backflow is installed in the main plant (as mentioned in Heymann and Rühling, 2016).

Conclusion

The present paper presented the operational behavior of a new prototype of a solar-fed bidirectional substation. This represents the final step of the THERMOSS project in which specifications, modeling, sizing, manufacturing and testing of such a prototype was performed. The prototype was sized for about 100kW in consumption and 70kW in production and can be considered in a first approach as a superposition of a consumption substation (with primary separation between DHW and SH) and a solar-fed feed-in substation.

The testing setup was composed of i) an emulated building SH and DHW consumption (real-time TRNSYS simulation associated to a Labview program and a set of thermo-hydraulic modules), ii) an emulated solar field (real-time TRNSYS simulation associated to a Labview program and a thermo-hydraulic module), and iii) a real 200m long linear DHN. The testing was performed during 12 days representative in terms of boundary conditions of the 12 months of the year. A set of indicators showed i) an inappropriate behavior of the controller during the coldest days, ii) various improvements of the solar fraction during the hottest days (with up to 52% of total solar fraction) and iii) rather good performances of the controller. The analysis was complemented by the observation of detailed dynamic result for one day of the sequence. Different operating regimes (full DHN consumption, solar self-consumption complemented by the DHN, solar self-consumption and feed-in) could be observed during a single day highlighting the modularity of the BSST to switch automatically (by hydraulic balancing) from one mode to the other. Unwanted low frequency oscillations between the modes were also observed.

As a conclusion, we can say that the experimental results shown here represent a satisfactory first step for such a bidirectional substation. However, there exists still some room for improvements regarding the control, especially by benefiting from the local consumer information to adapt the control of the feed-in loop. This is currently under investigation.

Fig. 1 :

 1 Fig. 1: Technical drawing of the BSST prototype

Fig. 2 :

 2 Fig. 2: Main operating principlesRegarding the substation control, the solar field pump (P_BSST_SOL) is controlled to obtain a given flow rate set point. Similarly to[START_REF] Lennermo | Control of decentralised solar district heating[END_REF], the latter is calculated using a cascaded combination of a theoretical solar field model (function of solar irradiation, return temperature TT_BSST_11 and ambient temperature) and a closed-loop PID leading to a desired set point for the outlet solar field temperature (TT_BSST_10). For the latter set point, an adapted heating law as a function of the ambient temperature is implemented. The feed-in pump (P_BSST_FEED) is controlled to have the same flow-rate as in the solar loop (accounting for the calorific value correction between the glycol on the solar loop and the water in the BSST loop). The extra set of valves (V_BSST_BP1 and V_BSST_BP2) are used for specific transient conditions. Based on the demand side and solar production, the flow rate will automatically be distributed, due to hydraulic balancing, such as in mode 1 or mode 2.

Fig. 3 :

 3 Fig. 3: Picture of the bidirectional substation with some components of Fig. 1 highlighted

Fig. 4 :

 4 Fig. 4: Dynamic whole system test methodology

Fig. 5 :

 5 Fig. 5: BSST test setup at CEA-INES: emulated building, emulated solar field and real DHN

4. 1 .Fig. 6 :

 16 Fig. 6: Schematic representation and nomenclature associated to thermal power fluxes (Q) around the substation Tab. 1: Energy (E) based indicators definition

  Feed-in energy divided by the summed fluxes at the entrance of the SST on the DHN side Substation Prosumer Rate with respect to DHN τ pro,DHN = E SST,feed-in E SST,DHN + E SST,feed-in Locally feed-in energy divided by summed fluxes at the entrance of all the SST of the entire DHN in,E,2K = ∑ P feed,in . ∆t| F BSST_DHN <0 & TT01>T sp,TT01 -2K ∑ P feed,in . ∆t| F BSST_DHN <0

Fig. 7 :Fig. 8 :

 78 Fig. 7: Experimental results obtained for day 6 of the sequence (Top: heat power; Bottom: Temperatures)

Fig. 9 :

 9 Fig. 9: Qualitative impact of the BSST on the real DHN pressure cone

Tab. 3: Indicators for the 12-days test sequence

	Day	𝛈 𝐒𝐅	𝛈 𝐒𝐏	𝛈 𝐭𝐡	𝐒𝐅	𝐒𝐅 𝐭𝐨𝐭	𝛕 𝐟𝐞𝐞𝐝	𝛕 𝐩𝐫𝐨	𝛕 𝐩𝐫𝐨,𝐃𝐇𝐍 𝛔 𝐟𝐞𝐞𝐝-𝐢𝐧 𝛆 𝐟𝐞𝐞𝐝-𝐢𝐧,𝐄,𝟐𝐊

Acknowledgments

The project THERMOSS has received funding from the European Union's Horizon 2020 research and innovation program under grant agreement No 723562.


Thomas Apel
email: thomas.apel@unibw.de

Serge Nicaise
email: snicaise@univ-valenciennes.fr

Johannes Pfefferer
email: johannes.pfefferer@unibw.de

Discretization of the Poisson equation with non-smooth data and emphasis on non-convex domains *

Keywords: Elliptic boundary value problem, method of transposition, very weak formulation, finite element method, discretization error estimate AMS subject classification 65N30; 65N15

Discretization of the Poisson equation with non-smooth data and emphasis on non-convex domains

Introduction

The motivation for this paper is to consider the boundary value problem -∆y = f in Ω, y = u on Γ := ∂Ω, (

with right hand side f ∈ H -1 (Ω) and boundary data u ∈ L 2 (Γ). We assume Ω ⊂ R 2 to be a bounded polygonal domain with boundary Γ. Such problems arise in optimal control when the Dirichlet boundary control is considered in L 2 (Γ) only, see for example the papers by Deckelnick, Günther, and Hinze, [START_REF] Deckelnick | Finite element approximation of Dirichlet boundary control for elliptic PDEs on two-and three-dimensional curved domains[END_REF], French and King, [START_REF] French | Approximation of an elliptic control problem by the finite element method[END_REF], May, Rannacher, and Vexler, [START_REF] May | Error analysis for a finite element approximation of elliptic Dirichlet boundary control problems[END_REF], and Apel, Mateos, Pfefferer, and Rösch, [START_REF] Apel | On the regularity of the solutions of Dirichlet optimal control problems in polygonal domains[END_REF]. On the continuous level we even admit more irregular data.

In Section 2 we analyze several ways how to understand the solution of the boundary value problem (1.1) for which we cannot expect a weak solution y ∈ H 1 (Ω). The most popular method to solve problem (1.1) is the transposition method that goes back to Lions and Magenes [START_REF] Lions | Problèmes aux limites non homogènes et applications[END_REF] and that is based on the use of some integration by parts. This formally leads to the very weak formulation: Find y ∈ Y such that (y, ∆v

) Ω = (u, ∂ n v) Γ -(f, v) Ω ∀v ∈ V (1.2)
with (w, v) G := G wv denoting the L 2 (G) scalar product or an appropriate duality product. The main issue is to find the appropriate trial space Y and test space V . In the convex case, it turns out that a good choice is Y = L 2 (Ω) and V = H 2 (Ω) ∩ H 1 0 (Ω). This case is investigated by many authors, some of them are mentioned in Subsection 2.1, but we will see that such a choice is not appropriate in the non-convex case (loss of uniqueness) and present some remedies (enlarged test spaces) in Subsections 2.2 and 2.3. We will further show in Subsection 2.2 that in the correct setting the very weak solution corresponds to the one obtained by an integral equation technique and hence has H 1/2 (Ω) regularity.

The main drawback of the very weak formulation is the fact that a conforming discretization of the test space should be made by C 1 -elements. Hence Berggren proposed in [START_REF] Berggren | Approximations of very weak solutions to boundary-value problems[END_REF] to introduce two new variables, ϕ := -∆v and ζ := ∂ n v allowing to perform a simpler numerical analysis, see Subsection 2.4. In Subsection 2.5 we propose another method that consists of regularizing the Dirichlet datum u by approximating it by a sequence of functions u h in H 1/2 (Γ) using for example an interpolation operator. This allows to compute a sequence of weak solutions y h ∈ H 1 (Ω), and we show that they converge to the very weak solution with an explicit convergence rate.

A negative result about the well-posedness of the weak formulation with L 1 (Γ)-data completes the discussion on the continuous level.

Section 3 is devoted to the numerical analysis. We start with Berggren's numerical approach and recall his error estimates in Subsection 3.1 for completeness. Next we perform in Subsection 3.2 a numerical analysis of our regularization approach and prove error estimates for the piecewise linear approximation on a family of conforming, quasi-uniform finite element meshes. Notice that it turns out that on the discrete level Berggren's approach is a particular case of our regularization strategy. The convergence order is 1 2 in the convex case but smaller in the non-convex case. This reduction can be explained by the singular behaviour of the solution of the dual problem. In our paper [START_REF] Apel | A dual singular complement method for the numerical solution of the Poisson equation with L 2 boundary data in non-convex domains[END_REF] we investigate the singular complement method to remedy the suboptimality of the standard finite element method in non-convex domains.

Finally, in Section 4 we present numerical tests in order to illustrate that our error estimates are sharp. The paper ends with some remarks about the three-dimensional case and about data with different regularity than assumed above.

Analysis of the boundary value problem

In this section we analyze several ways how to understand the solution of the boundary value problem (1.1) for which we cannot expect a weak solution y ∈ H 1 (Ω).

For keeping the notation succinct we assume that the polygonal domain Ω has at most one non-convex corner with interior angle, denoted by ω. Let r, θ be the corresponding polar coordinates and define λ := π/ω. The boundary segments of Ω are denoted by Γ j , j = 1, . . . , N .

Method of transposition in convex polygonal domains

The method of transposition goes back at least to Lions and Magenes, [START_REF] Lions | Problèmes aux limites non homogènes et applications[END_REF], and is used by several other authors including French and King, [START_REF] French | Approximation of an elliptic control problem by the finite element method[END_REF], Casas and Raymond, [START_REF] Casas | Error estimates for the numerical approximation of Dirichlet boundary control for semilinear elliptic equations[END_REF], Deckelnick, Günther, and Hinze, [START_REF] Deckelnick | Finite element approximation of Dirichlet boundary control for elliptic PDEs on two-and three-dimensional curved domains[END_REF], and May, Rannacher, and Vexler, [START_REF] May | Error analysis for a finite element approximation of elliptic Dirichlet boundary control problems[END_REF]. Since by partial integration the derivation

(f, v) Ω = -(∆y, v) Ω = (∇y, ∇v) Ω for v ∈ H 1 0 (Ω) = (y, ∂ n v) Γ -(y, ∆v) Ω for v ∈ H 2 (Ω) ∩ H 1 0 (Ω)
is valid we get the very weak formulation: Find

y ∈ L 2 (Ω) : (y, ∆v) Ω = (u, ∂ n v) Γ -(f, v) Ω ∀v ∈ H 2 (Ω) ∩ H 1 0 (Ω). (2.1)
May, Rannacher, and Vexler, [START_REF] May | Error analysis for a finite element approximation of elliptic Dirichlet boundary control problems[END_REF], proved the existence of a solution in L 2 (Ω) in the case of a convex polygonal domain Ω.

Lemma 2.1. If the domain Ω is convex there exists a unique solution y ∈ L 2 (Ω) of problem (2.1) that satisfies the a priori error estimate

y L 2 (Ω) ≤ c u N j=1 H 1/2 00 (Γ j ) + f (H 2 (Ω)∩H 1 0 (Ω)) (2.2) provided that u ∈ N j=1 H 1/2 00 (Γ j ) and f ∈ H 2 (Ω) ∩ H 1 0 (Ω) .
Recall that H 1/2 00 (Γ j ) is the space of functions whose extension by zero to Γ is in

H 1/2 (Γ).
Proof. For being self-contained we sketch here the proof of [START_REF] May | Error analysis for a finite element approximation of elliptic Dirichlet boundary control problems[END_REF]Lemma 2.1]. The idea is first to assume more regular data, u ∈ H 1/2 (Γ), f ∈ H -1 (Ω) such that a weak solution y ∈ H 1 (Ω) exists, then to show (2.2), and finally to use a standard density argument since

H 2 (Ω) ∩ H 1 0 (Ω) c → H 1 0 (Ω) and H 1/2 (Γ) c → L 2 (Γ) = N j=1 L 2 (Γ j ) c → N j=1 H 1/2 00 (Γ j ) . The estimate (2.
2) is proven by using a duality argument. Due to the convexity of the domain there exists a solution w ∈ H 2 (Ω) of the auxiliary problem

-∆w = y in Ω, w = 0 on Γ, (2.3)
such that

y 2 L 2 (Ω) = (y, -∆w) Ω = (f, w) Ω -(u, ∂ n w) Γ ≤ c u N j=1 H 1/2 00 (Γ j ) + f (H 2 (Ω)∩H 1 0 (Ω)) w H 2 (Ω) ,
where one uses that the mapping [START_REF] Grisvard | Elliptic problems in nonsmooth domains[END_REF]Thm. 1.5.2.8]. The desired estimated (2.2) is then obtained by using the a priori estimate w H 2 (Ω) ≤ c y L 2 (Ω) and by division by y L 2 (Ω) .

H 2 (Ω) ∩ H 1 0 (Ω) → N j=1 H 1/2 00 (Γ j ), w → ∂ n w, is surjective due to
The method of proof of this lemma will be revisited in Section 2.5, where we start the discussion of the regularization approach. However, there we will work with different function spaces such that non-convex domains are included in the theory as well. The given proof of Lemma 2.1 is even restricted to convex domains since the isomorphism

∆w ∈ L 2 (Ω), w |Γ = 0 ⇔ w ∈ H 2 (Ω) ∩ H 1 0 (Ω)
is used. If the domain is non-convex one loses at least the uniqueness of the solution y of (2.1). For example, take Ω = {(r cos θ, r sin θ) ∈ R 2 : 0 < r < 1, 0 < θ < ω} and λ = π/ω, then both y 1 = r λ sin(λθ) and y 2 = r -λ sin(λθ) are harmonic in Ω and y 1 = y 2 on Γ. Both satisfy (2.1) with f ≡ 0 and g = y 1 = y 2 on Γ. Hence, one needs a larger test space in order to rule out y 2 .

Method of transposition in general polygonal domains

In a first instance we replace the test space H 2 (Ω) ∩ H 1 0 (Ω) by

V := H 1 0 (Ω) ∩ H 1 ∆ (Ω) with H 1 ∆ (Ω) := {v ∈ H 1 (Ω) : ∆v ∈ L 2 (Ω)}. (2.4) Since |v| H 1 (Ω) ≤ c ∆v L 2 (Ω) , the graph norm in V , that is ∆v L 2 (Ω) + |v| H 1 (Ω) , is equivalent to ∆v L 2 (Ω)
such that we will use henceforth

v V = ∆v L 2 (Ω) .
Furthermore, let us denote by V Γ the space of normal derivatives ∂ n v of functions v ∈ V . According to [START_REF] Grisvard | Elliptic problems in nonsmooth domains[END_REF]Theorem 1.5.3.10] this space is well-defined and a subspace of N j=1 H

1/2 00 (Γ j ) . The natural norm in V Γ is given by

g V Γ := inf { v V : v ∈ V, ∂ n v = g} .
In particular, we have for v ∈ V that

∂ n v V Γ ≤ v V . (2.5)
Since the previous definitions of the spaces V and V Γ are rather formal let us discuss the structure of these spaces.

Remark 2.2. The spaces V and V Γ can be characterized as follows:

1. If Ω is convex the spaces V and H 2 (Ω) ∩ H 1 0 (Ω) coincide. However, in non-convex domains the situation is different. In this case there is the splitting

V = H 2 (Ω) ∩ H 1 0 (Ω) ⊕ Span{ξ(r)
V Γ =   N j=1 H 1/2 00 (Γ j )   ⊕ Span{ξ(r) r λ-1 }. 3. In the non-convex case there is V → H s (Ω) ∩ H 1 0 (Ω) for s < 1 + λ and V Γ → H t (Γ) for t < λ -1 2 . This implies H s (Ω) ∩ H 1 0 (Ω) → V and H t (Γ) → V Γ . Lemma 2.3. Let f ∈ V and u ∈ V Γ .
Then there exists a unique solution

y ∈ L 2 (Ω) : (y, ∆v) Ω = (u, ∂ n v) Γ -(f, v) Ω ∀v ∈ V (2.6) with y L 2 (Ω) ≤ u V Γ + f V .
Proof. The proof of this lemma is based on the Babuška-Lax-Milgram theorem. Due to (2.5) the right hand side of (2.6) defines a linear functional on V . Moreover, the bilinear form is bounded on L 2 (Ω) × V . The inf-sup conditions are proved by using the isomorphism ∆v

∈ L 2 (Ω), v |Γ = 0 ⇔ v ∈ V.
In particular, we obtain by taking y = ∆v sup

y∈L 2 (Ω) |(y, ∆v) Ω | y L 2 (Ω) ≥ (∆v, ∆v) Ω ∆v L 2 (Ω) = ∆v L 2 (Ω) = v V ,
and by taking the solution v ∈ V of ∆v = y

sup v∈V |(y, ∆v) Ω | v V ≥ (y, y) Ω y L 2 (Ω) = y L 2 (Ω) . (2.7)
The existence of the unique solution y ∈ L 2 (Ω) of problem (2.6) follows from the standard Babuška-Lax-Milgram theorem, see for example [START_REF] Babuška | Error-bounds for finite element method[END_REF]Theorem 2.1]. The a priori estimate follows from (2.7), and (2.6) and (2.5),

y L 2 (Ω) ≤ sup v∈V |(y, ∆v) Ω | v V = sup v∈V |(u, ∂ n v) Γ -(f, v) Ω | v V ≤ u V Γ + f V .
Note that if a weak solution y ∈ {v ∈ H 1 (Ω) : v| Γ = u} exists, then with the help of the Green formula

(∂ n v, χ) Γ = (∇v, ∇χ) Ω + (∆v, χ) Ω ∀v ∈ H 1 ∆ (Ω) ⊃ V, ∀χ ∈ H 1 (Ω), (2.8)
see Lemma 3.4 in the paper [START_REF] Costabel | Boundary integral operators on Lipschitz domains: elementary results[END_REF] by Costabel, it is also a very weak solution. (Set χ = y and use (∇v, ∇y) Ω = (f, v) Ω .) A possible difficulty with this formulation is that a conforming discretization with a finite-dimensional space V h ⊂ V would require the use of C 1 -functions. This is simple for one-dimensional domains Ω but requires a lot of degrees of freedom in two (and more) dimensions.

We finish this subsection with a regularity result.

Lemma 2.4. The unique solution y ∈ L 2 (Ω) of problem (2.6) with u ∈ L 2 (Γ) and f = 0 belongs to H 1/2 (Ω) and to

W 1,2 (Ω) := {z ∈ L 2 (Ω) : δ 1/2 ∇z ∈ L 2 (Ω) 2 } (2.9)
where δ(x) is the distance of x to the boundary Γ. Furthermore, there exists a positive constant c such that

y H 1/2 (Ω) + δ 1/2 ∇y L 2 (Ω) 2 ≤ c u L 2 (Γ) .
(2.10)

Remark 2.5. The book by Chabrowski, [START_REF] Chabrowski | The Dirichlet problem with L 2 -boundary data for elliptic linear equations[END_REF], deals exclusively with the Dirichlet problem with L 2 boundary data for elliptic linear equations. The solution is searched there in the Sobolev space (2.9) but domains of class C 1,1 were considered only.

Proof of Lemma 2.4. In a first step, we use an integral representation and some properties of the layer potentials to get a solution with the appropriate regularity. Theorem 4.2 of Verchota's paper [START_REF] Verchota | Layer potentials and regularity for the Dirichlet problem for Laplace's equation in Lipschitz domains[END_REF] shows that the operator 1 2 I + K, K being the boundary double layer potential, is an isomorphism from L 2 (Γ) into itself (see also Corollary 4.5 of [START_REF] Mitrea | Boundary layer methods for Lipschitz domains in Riemannian manifolds[END_REF]). According to the trace property for the double layer potential K (see Section 1 and Corollary 3.2 of [START_REF] Verchota | Layer potentials and regularity for the Dirichlet problem for Laplace's equation in Lipschitz domains[END_REF]), there exists a unique harmonic function z such that z → u a.e. in nontangential cones with the notation from [START_REF] Verchota | Layer potentials and regularity for the Dirichlet problem for Laplace's equation in Lipschitz domains[END_REF], and is given by

z = K( 1 2 I + K) -1 u. (2.11)
But due to the above mentioned isomorphism property and Theorem 1 of Costabel's paper [START_REF] Costabel | Boundary integral operators on Lipschitz domains: elementary results[END_REF] (and its following Remark), we obtain that

z H 1/2 (Ω) ≤ c u L 2 (Γ) .
(2.12)

Note that Theorems 5.3, 5.4 and Corollary 5.5 of the paper [START_REF] Jerison | The inhomogeneous Dirichlet problem in Lipschitz domains[END_REF] by Jerison and Kenig also yield

δ 1/2 ∇z L 2 (Ω) 2 ≤ c u L 2 (Γ) . (2.13)
The second step is to show that z is the very weak solution, hence by uniqueness, we will get y = z. The regularity of y and the estimate (2.10) follow from our first step. For that last purpose, we use a density argument. Indeed let u n ∈ H 1 (Γ) be a sequence of functions such that

u n → u in L 2 (Γ), as n → ∞. (2.14) Consider z n = K( 1 2 I + K) -1
u n and let y n ∈ L 2 (Ω) be the unique solution of (2.6) with boundary datum u n and right hand side f = 0 (that is in H 1 (Ω)). Then by the estimate (2.12) and Lemma 2.3, we get

y n → y in L 2 (Ω), as n → ∞, z n → z in H 1/2 (Ω), as n → ∞.
Furthermore by Theorem 5.15 of [START_REF] Jerison | The inhomogeneous Dirichlet problem in Lipschitz domains[END_REF] z n satisfies

γz n = u n on Γ,
where γ is the trace operator from H 1 (Ω) into H 1/2 (Γ). Hence we directly deduce that y n = z n and by the above convergence property we conclude that y = z.

Corollary 2.6. The unique solution y ∈ L 2 (Ω) of problem (2.6) with u ∈ L 2 (Γ) and f ∈ H -1 (Ω) belongs to H 1/2 (Ω) and to W 1,2 (Ω) from (2.9). There exists a positive constant c such that

y H 1/2 (Ω) + δ 1/2 ∇y L 2 (Ω) 2 ≤ c u L 2 (Γ) + f H -1 (Ω) .

Method of transposition employing weighted Sobolev spaces

Alternatively to the space

H 1 ∆ (Ω) ∩ H 1 0 (Ω), one can use the test space V 2,2 β (Ω) ∩ H 1 0 (Ω), β > 1 - π ω ,
in the non-convex case, where V 2,2 β (Ω) is a weighted Sobolev space of the class

V k,p β (Ω) := v ∈ D (Ω) : v V k,p β (Ω) < ∞ , v p V k,p β (Ω) := |α|≤k Ω |r β-k+|α| D α v| p , (2.15)
where we use standard multi-index notation. For later use we also introduce [START_REF] Nicaise | Polygonal Interface Problems, volume 39 of Methoden und Verfahren der mathematischen Physik[END_REF]Theorem 1.31]. In the next lemma, we will use the spaces

L 2 β (Ω) := V 0,2 β (Ω). First derivatives of V 2,2 β (Ω)-functions belong to V 1,2 β (Ω) by definition. The trace space of V 1,2 β (Ω) is N j=1 V 1/2,2 β (Γ j ), see [18, Lemma 1.2] or
V β := V 2,2 β (Ω) ∩ H 1 0 (Ω) for ω > π, H 2 (Ω) ∩ H 1 0 (Ω) for ω < π, Y β := L 2 -β (Ω) for ω > π, L 2 (Ω) for ω < π, for β ∈ (1 -π ω , 1]
. We endow V β with the V 2,2 β (Ω)-norm for ω > π and the H 2 (Ω)-norm for ω < π, as well as Y β with the L 2 -β (Ω)-norm for ω > π and the L 2 (Ω)-norm for ω < π.

Remark 2.7. Let us discuss the definition of the space V β and the restriction of the weight β to the interval (1

-π ω , 1] in the non-convex case: 1. We require β ∈ (1 -π ω , 1] in order to have the isomorphism ∆v ∈ L 2 β (Ω), v |Γ = 0 ⇔ v ∈ V β . (2.16)
2. It is possible to use the test space V 2,2 β (Ω) ∩ H 1 0 (Ω) in convex domains as well. However, this implies a loss of information about the solution since this test space is smaller than H 2 (Ω) ∩ H 1 0 (Ω) due to the fact that the weight β can be negative.

3. Note that V 1/2,2 β (Γ j ) → L 2 (Γ j ) for β ≤ 1 2 . This implies L 2 (Γ) = N j=1 L 2 (Γ j ) → N j=1 V 1/2,2 β (Γ j )
. This means that L 2 -boundary data are included in the fol-

lowing discussion if β ≤ 1 2 . Lemma 2.8. Let β ∈ (1 -π ω , 1] and assume that f ∈ V β and u ∈ N j=1 V 1/2,2 β (Γ j ) .
Then there exists a unique solution

y ∈ Y β → L 2 (Ω) : (y, ∆v) Ω = (u, ∂ n v) Γ -(f, v) Ω ∀v ∈ V β (2.17) with y Y β ≤ c f V β + u N j=1 V 1/2,2 β (Γ j )
.

Proof. The convex case was already treated in Lemma 2.1, hence we focus on the nonconvex case. We proceed as in Lemma 2.3. The right hand side of (2.17) defines a continuous functional on V β . Furthermore, the bilinear form is bounded on

Y β × V β , |(y, ∆v) Ω | ≤ r -β y L 2 (Ω) r β ∆v L 2 (Ω) ≤ y Y β v V β .
The inf-sup conditions for the bilinear form are proved by using the isomorphism (2.16), in particular v V β ≤ c ∆v L 2 β (Ω) ; we obtain by taking y = r 2β ∆v sup

y∈Y β |(y, ∆v) Ω | y Y β ≥ (r β ∆v, r β ∆v) Ω r 2β ∆v L 2 -β (Ω) = ∆v 2 L 2 β (Ω) ∆v L 2 β (Ω) ≥ c -1 v V β ,
and by taking the solution

v ∈ V β of ∆v = r -2β y with v V β ≤ c ∆v L 2 β (Ω) = c y L 2 -β (Ω) sup v∈V β |(y, ∆v) Ω | v V β ≥ (r -β y, r -β y) Ω c y L 2 -β (Ω) = c -1 y L 2 -β (Ω) . (2.18)
The existence of the unique solution y ∈ L 2 -β (Ω) of problem (2.17) follows now from the standard Babuška-Lax-Milgram theorem, see for example [START_REF] Babuška | Error-bounds for finite element method[END_REF]Theorem 2.1]. The a priori estimate is obtained with (2.18) and (2.17),

y Y β ≤ c sup v∈V β |(y, ∆v) Ω | v V β = c sup v∈V β |(u, ∂ n v) Γ -(f, v) Ω | v V β .
This ends the proof since we already noticed that the enumerator defines a continuous functional on V β .

Berggren's approach

Berggren's approach [START_REF] Berggren | Approximations of very weak solutions to boundary-value problems[END_REF] avoids test functions in H 1 ∆ (Ω) ∩ H 1 0 (Ω) in an explicit way. It can be explained as if we substitute ϕ := -∆v and ζ := ∂ n v in (2.6),

y ∈ L 2 (Ω) : (y, ϕ) Ω = -(u, ζ) Γ + (f, v) Ω ∀ϕ ∈ L 2 (Ω).
(2. [START_REF] Mitrea | Boundary layer methods for Lipschitz domains in Riemannian manifolds[END_REF] The relationship between ϕ ∈ L 2 (Ω) and both v ∈ V and ζ ∈ V Γ can be expressed by the weak formulation of the Poisson equation,

v ∈ V : (∇v, ∇ψ) Ω = (ϕ, ψ) Ω ∀ψ ∈ H 1 0 (Ω) (2.20)
and a reformulation of the Green formula (2.8) in the form

ζ ∈ V Γ : (ζ, χ) Γ = (∇v, ∇χ) Ω -(ϕ, χ) Ω ∀χ ∈ H 1 (Ω) \ H 1 0 (Ω). ( 2

.21)

Note that Berggren's formulation is not a system with three unknown functions since the second and third equations compute actions on the test function ϕ. Indeed, let S : L 2 (Ω) → V and F : L 2 (Ω) → V Γ be the solution operators of (2.20) and (2.21), respectively, defined by Sϕ := v and F ϕ := ζ, then we could also write

y ∈ L 2 (Ω) : (y, ϕ) Ω = -(u, F ϕ) Γ -(f, Sϕ) Ω ∀ϕ ∈ L 2 (Ω) instead of (2.
-(y, ∆v) Ω = -(u, ∂ n v) Γ + (f, v) Ω ∀v ∈ V
due to the isometry between L 2 (Ω) and V .

Remark 2.10. Berggren used the regularity v ∈ H 3/2+ (Ω) with some > 0 which implies ζ ∈ H (Γ). However, he did not consider the maximal domain of the elliptic operator, i.e., v ∈ V and ζ ∈ V Γ . But with the explanations of Subsection 2.2 these regularities should be obvious. Thus the result of Lemma 2.9 is slightly more general than that of Berggren.

The regularization approach

A further idea is to regularize the boundary data and then to apply standard methods. This approach has already been considered within the proof of Lemma 2.1. In contrast, we do not use the isomorphism

∆w ∈ L 2 (Ω), w |Γ = 0 ⇔ w ∈ H 2 (Ω) ∩ H 1 0 (Ω),
which can only be employed in case of convex domains, but the isomorphism

∆v ∈ L 2 (Ω), v |Γ = 0 ⇔ v ∈ V.
This allows us to apply the regularization approach in the non-convex case as well. Moreover, we propose two different strategies how the regularized Dirichlet boundary data can be constructed in an explicit way. Thereby we will be able in Subsection 3.2 to calculate approximate solutions of the regularized problems based on a finite element method. For the data we assume henceforth u ∈ L 2 (Γ) and f ∈ H -1 (Ω). This is not only for simplicity but also due to the fact that already for Dirichlet boundary data in L 2 (Γ) the convergence rates of the approximate solutions in Subsection 3.2 tend to zero as the maximal interior angle tends to 2π. We start with general convergence results for the regularized solutions. To this end let u h ∈ H 1/2 (Γ) be a sequence of functions such that

lim h→0 u -u h L 2 (Γ) = 0. Let now y h ∈ Y h * := {v ∈ H 1 (Ω) : v| Γ = u h } be the variational solution, y h ∈ Y h * : (∇y h , ∇v) Ω = (f, v) Ω ∀v ∈ H 1 0 (Ω). (2.22)
Lemma 2.11. Let u ∈ L 2 (Γ) and f ∈ H -1 (Ω). Then the limit y := lim h→0 y h exists, belongs to L 2 (Ω), and is the very weak solution, that means it satisfies (2.6).

Proof. First we show that y h is a Cauchy sequence in L 2 (Ω). From (2.22) and Green's formula, we have for any v ∈ V ,

(f, v) Ω = (∇y h , ∇v) Ω = -(y h , ∆v) Ω + (y h , ∂ n v) Γ , (f, v) Ω = (∇y h , ∇v) Ω = -(y h , ∆v) Ω + (y h , ∂ n v) Γ .
Hence due to y h = u h and y h = u h on Γ, we deduce that

(y h -y h , ∆v) Ω = (u h -u h , ∂ n v) Γ ∀v ∈ V. (2.23) Now for any z ∈ L 2 (Ω), let v z ∈ V be such that ∆v z = z, (2.24)
that clearly satisfies

∂ n v z L 2 (Γ) ≤ c v z H s (Ω) ≤ c z L 2 (Ω) (2.25)
with some s ∈ 3 2 , 1 + λ , s ≤ 2. Finally, we obtain with (2.23) and (2.25)

y h -y h L 2 (Ω) = sup z∈L 2 (Ω),z =0 (y h -y h , z) Ω z L 2 (Ω) = sup z∈L 2 (Ω),z =0 (u h -u h , ∂ n v z ) Γ z L 2 (Ω) ≤ u h -u h L 2 (Γ) sup z∈L 2 (Ω),z =0 ∂ n v z L 2 (Γ) z L 2 (Ω) = c u h -u h L 2 (Γ) .
Since u h converges in L 2 (Γ), it is a Cauchy sequence and hence also y h is a Cauchy sequence and converges in L 2 (Ω) by the completeness of L 2 (Ω). From V ⊂ H 1 0 (Ω) we obtain by (2.22) and the Green formula (2.8)

(f, v) Ω = (∇y h , ∇v) Ω = (∆v, y h ) Ω -(∂ n v, y h ) Γ = (∆v, y h ) Ω -(∂ n v, u h ) Γ ∀v ∈ V.
Since ∆v ∈ L 2 (Ω) and ∂ n v ∈ L 2 (Γ) we can pass to the limit and obtain that the limit function y satisfies (2.6).

We can estimate the regularization error by a similar technique.

Lemma 2.12.

Let s = 1 2 if Ω is convex and s ∈ [0, λ -1 2 ) if Ω is non-convex. Then the estimate y -y h L 2 (Ω) ≤ c u -u h H -s (Γ)
holds.

Proof. We use the approach of the proof of Lemma 2.11 and just replace u h by u and y h by y to get

(y -y h , ∆v) Ω = (u -u h , ∂ n v) Γ ∀v ∈ V. (2.26)
Again, for any z ∈ L 2 (Ω), we let v z ∈ V be such that ∆v z = z but estimate now in a sharper way

∂ n v z H s (Γ) ≤ c v z H s+3/2 (Ω) ≤ c z L 2 (Ω) (2.27)
As in the previous proof we get

y -y h L 2 (Ω) = sup z∈L 2 (Ω),z =0 (u -u h , ∂ n v z ) Ω z L 2 (Ω) ≤ u -u h H -s (Γ) sup z∈L 2 (Ω),z =0 ∂ n v z H s (Γ) z L 2 (Ω) ≤ c u -u h H -s (Γ) .
Actually, the proof is for s > 0 but of course the statement holds when s is decreased.

A choice for the construction of the regularized function u h could be the use of the L 2 (Γ)-projection Π h u into a piecewise polynomial space on the boundary (which we call Y ∂ h in Section 3) or the use of the Carstensen interpolant C h u, see [START_REF] Carstensen | Quasi-interpolation and a posteriori error analysis in finite element methods[END_REF]. Namely, if N Γ is the set of nodes of the triangulation on the boundary, we set

C h u = x∈N Γ π x (u)λ x ,
where λ x is the standard hat function related to x and

π x (u) = Γ uλ x Γ λ x = (u, λ x ) Γ (1, λ x ) Γ .
The advantages of the interpolant in comparison with the L 2 -projection are its local definition and the property

u ∈ [a, b] ⇒ C h u ∈ [a, b],
see [START_REF] De Los Reyes | Finite element error analysis for stateconstrained optimal control of the Stokes equations[END_REF]; a disadvantage for our application in optimal control is that C h u h = u h for piecewise linear u h . We prove now regularization error estimates for the case that the regularized function u h is constructed via C h u or Π h u.

Lemma 2.13. If u h is the piecewise linear Carstensen interpolant of u or the L 2 (Γ)projection of u into a space of piecewise linear functions, then there holds

u -u h H -s (Γ) ≤ ch s u L 2 (Γ) , s ∈ [0, 1]
as well as

u h L 2 (Γ) ≤ c u L 2 (Γ) .
Proof. The interpolation error estimate and stability result are derived in [START_REF] Carstensen | Quasi-interpolation and a posteriori error analysis in finite element methods[END_REF][START_REF] De Los Reyes | Finite element error analysis for stateconstrained optimal control of the Stokes equations[END_REF] for domains. The proofs can be transferred to estimates on the boundary Γ. For the sake of completeness we sketch these proofs in the Appendix, see Lemma A.2.

In the case of the L 2 -projection, the second estimate holds with constant one. For the first estimate, we notice that by using the properties of the L 2 -projection, we have

(u -Π h u, ϕ) Γ = (u -Π h u, ϕ -Π h ϕ) Γ = (u, ϕ -Π h ϕ) Γ ≤ u L 2 (Γ) ϕ -Π h ϕ L 2 (Γ) ≤ u L 2 (Γ) ϕ -C h ϕ L 2 (Γ) ≤ ch s u L 2 (Γ) ϕ H s (Γ)
where we used Lemma A.2 in the last step. We conclude

u -u h H -s (Γ) = sup ϕ∈H s (Γ),ϕ =0 (u -Π h u, ϕ) Γ ϕ H s (Γ) ≤ ch s u L 2 (Γ) .
which is the assertion.

By setting s = 0 in the previous lemma we obtain

u -u h L 2 (Γ) ≤ c u L 2 (Γ)
for the different choices of the regularized function u h . This means that for u ∈ L 2 (Γ) the difference u -u h is uniformly bounded in L 2 (Γ) independent of h. However, we require strong convergence in

L 2 (Γ) for u ∈ L 2 (Γ), i.e. lim h→0 u -u h L 2 (Γ) = 0.
This is subject of the next lemma. A comparable result for the Ritz-projection can be found in e.g. [8, Theorem 3.2.3].

Lemma 2.14. Let u ∈ L 2 (Γ) and let u h be the piecewise linear Carstensen interpolant of u or the L 2 (Γ)-projection of u into a space of piecewise linear functions. Then there holds lim

h→0 u -u h L 2 (Γ) = 0.
Proof. We show the validity of this lemma for the Carstensen interpolant. The convergence result for the L 2 (Γ)-projection can be proven analogously. Due to the compact embedding H 1 (Γ) c → L 2 (Γ) there exists a sequence of functions

u n ∈ H 1 (Γ) such that lim n→∞ u -u n L 2 (Γ) = 0 with u n H 1 (Γ) ≤ c n , (2.28)
where the constant c n may depend on n. Thus, for every ε > 0 there is a positive integer

N such that u -u n L 2 (Γ) ≤ ε 2 (2.29)
for n ≥ N . By inserting the function u n and its Carstensen interpolant C h u n as intermediate functions into the desired term we obtain

u -C h u L 2 (Γ) ≤ u -u n L 2 (Γ) + u n -C h u n L 2 (Γ) + C h (u n -u) L 2 (Γ) .
The results of Lemma 2.13 and the inequalities (2.28) and (2.29) imply

u -C h u L 2 (Γ) ≤ c ε 2 + hc n .
Since for every constant c n there is a parameter h n such that hc n ≤ ε 2 for all h ≤ h n we arrive at

u -C h u L 2 (Γ) ≤ cε
and the desired result follows.

A negative result

Since the boundary datum u ∈ L 2 (Γ) → L 1 (Γ) is trace of a function w ∈ W 1,1 (Ω) an idea could be to search y = y 0 + w with y 0 | Γ = 0 and

y 0 ∈ W 1,1 0 (Ω) : (∇y 0 , ∇v) Ω = (f, v) Ω -(∇w, ∇v) Ω =: (F, v) Ω ∀v ∈ W 1,∞ 0 (Ω).
However, bilinear forms a(u, v) : W 1,p 0 (Ω) × W 1,p 0 (Ω) → C are investigated in [START_REF] Simader | On Dirichlet's boundary value problem. An L p -theory based on a generalization of G ȧrding's inequality[END_REF][START_REF] Simader | The Dirichlet problem for the Laplacian in bounded and unbounded domains[END_REF] and it is shown that

∀F ∈ W 1,p (Ω) * ∃z ∈ W 1,p 0 (Ω) : a(z, v) = (F, v) Ω ∀v ∈ W 1,p 0 (Ω) holds for p ∈ (1, ∞) only since c p → ∞ for p → 1 in the inf-sup condition c p sup φ∈W 1,p 0 (Ω) (∇z, ∇φ) Ω ∇φ L p (Ω) ≥ ∇z L p (Ω) for 1 < p < ∞.

Discretization of the boundary value problem

Let (T h ) h>0 be a family of conforming, quasi-uniform finite element meshes, and introduce the finite element spaces

Y h = {v h ∈ H 1 (Ω) : v h | T ∈ P 1 ∀T ∈ T h }, Y 0h = Y h ∩ H 1 0 (Ω), Y ∂ h = Y h | ∂Ω . (3.1)

Berggren's approach

Let u ∈ L 2 (Γ) and f ∈ H -1 (Ω). We discretize the formulation (2.19), (2.20), (2.21) in a straightforward manner,

y h ∈ Y h : (y h , ϕ h ) Ω = -(u, ζ h ) Γ + (f, v h ) Ω ∀ϕ h ∈ Y h , v h ∈ Y 0h : (∇v h , ∇ψ h ) Ω = (ϕ h , ψ h ) Ω ∀ψ h ∈ Y 0h , ζ h ∈ Y ∂ h : (ζ h , χ h ) Γ = (∇v h , ∇χ h ) Ω -(ϕ h , χ h ) Ω ∀χ h ∈ Y h \ Y 0h .
Note that ζ h = F ϕ h and v h = Sϕ h with F and S from Subsection 2.4, hence it is only an approximate Galerkin formulation. Berggren showed in Theorem 5.2 of [START_REF] Berggren | Approximations of very weak solutions to boundary-value problems[END_REF] that this formulation is equivalent with the standard finite element approximation with L 2 -projection of the boundary data,

y h ∈ Y h : (∇y h , ∇ϕ h ) Ω = (f, ϕ h ) Ω ∀ϕ h ∈ Y 0h , (3.2)
(y h , ϕ h ) Γ = (u, ϕ h ) Γ ∀ϕ h ∈ Y ∂ h . (3.3)
Note that y h = u h on Γ when u = u h ∈ Y ∂ h . This will be of interest in the discretization of optimal control problems.

Berggren proved also the following discretization error estimate, see Theorem 5.5 of [START_REF] Berggren | Approximations of very weak solutions to boundary-value problems[END_REF].

Lemma 3.1. Let ω be the maximal interior angle of the domain Ω, and denote by λ = π/ω the corresponding singularity exponent. Let s ∈ (0, 1 2 ] be a real number with s < λ - 1 2 , and let s ∈ [0, s ) be a further real number. Then the error estimate

y -y h L 2 (Ω) ≤ c h s y H s (Ω) + h s u L 2 (Γ) + h s +1/2 f H -1 (Ω) = O(h s )
holds, this means that we have convergence order s,

s = min 1 2 , λ -1 2 -ε = 1 2 -ε for convex domains, λ -1 2 -ε for non-convex domains, ε > 0 arbitrary.
Note that s → 0 for ω → 2π. We will show in the next section that in the convex case the convergence order is 1 2 , without ε.

Regularization approach

We consider a regularization strategy such that u h ∈ Y ∂ h , see Subsection 2.5. Recall that the corresponding solution

y h ∈ Y h * := {v ∈ H 1 (Ω) : v| Γ = u h } is defined via (2.22
). For a regularization using the L 2 (Γ)-projection or the Carstensen interpolant we have the regularization error estimate

y -y h L 2 (Ω) ≤ c u -u h H -s (Γ) ≤ ch s u L 2 (Γ) (3.4) with s = 1 2 if Ω is convex and s ∈ [0, λ -1 2 ) if Ω is non-convex,
y h ∈ Y * h : (∇y h , ∇v h ) Ω = (f, v h ) Ω ∀v h ∈ Y 0h . (3.5)
Note that, if we construct u h by the L 2 -projection, we recover the Berggren approach as a special case.

By the Poincaré inequality we obtain

∇y h 0 L 2 (Ω) ≤ c f H -1 (Ω) + ∇(B h u h ) L 2 (Ω) , ∇y 0h L 2 (Ω) ≤ c f H -1 (Ω) + ∇(B h u h ) L 2 (Ω) .
With the help of (3.6) we arrive at

∇y h 0 L 2 (Ω) + ∇y 0h L 2 (Ω) ≤ c f H -1 (Ω) + h -1/2 u L 2 (Γ) . (3.10)
Now as before we start with

y h 0 -y 0h L 2 (Ω) = sup z∈L 2 (Ω),z =0 (y h 0 -y 0h , z) Ω z L 2 (Ω) .
Letting again v z ∈ V be such that ∆v z = z, we get

y h 0 -y 0h L 2 (Ω) = sup z∈L 2 (Ω),z =0 (∇(y h 0 -y 0h ), ∇v z ) Ω z L 2 (Ω) ,
and therefore thanks to (3.9) we arrive at

y h 0 -y 0h L 2 (Ω) = sup z∈L 2 (Ω),z =0 (∇(y h 0 -y 0h ), ∇(v z -I h v z )) Ω z L 2 (Ω) ,
where I h is the standard Lagrange interpolation operator. By the Cauchy-Schwarz inequality, the well-known estimate

v z -I h v z 1,Ω ≤ ch 1/2+s v z H 3/2+s (Ω) ≤ ch 1/2+s z L 2 (Ω) ,
with s = 1 2 in the convex case and s < λ -1 2 in the non-convex case, and the estimate (2.27) we obtain y h 0 -y 0h L 2 (Ω) ≤ ch 1/2+s ∇(y h 0 -y 0h ) L 2 (Ω) . Using the a priori estimate (3.10), we conclude that

y h 0 -y 0h L 2 (Ω) ≤ ch s h 1/2 f H -1 (Ω) + u L 2 (Γ) .
Note that we have now proved the convergence order 1 2 for Berggren's approach, in the convex case. This section is devoted to the numerical verification of our theoretical results. For that purpose we present an example with known solution. Furthermore, to examine the influence of the corner singularities, we consider different polygonal domains Ω ω depending on an interior angle ω ∈ (0, 2π). These domains are defined by

Ω ω := (-1, 1) 2 ∩ {x ∈ R 2 : (r(x), ϕ(x)) ∈ (0, √ 2] × [0, ω]},
where r and ϕ stand for the polar coordinates located at the origin. The boundary of Ω ω is denoted by Γ ω which is decomposed into straight line segments Γ j , j = 1, . . . , m(ω), counting counterclockwise beginning at the origin. As numerical example we consider the problem -∆y = 0 in Ω ω , y = u on Γ j , j = 1, . . . , m(ω).

The boundary datum u is chosen as u := r -0.4999 sin(-0.4999ϕ) on Γ ω .

This function belongs to L p (Γ) for every p < 2.0004. The exact solution of our problem is simply y = r -0.4999 sin(-0.4999ϕ), since y is harmonic. A plot of it can be seen in Figure 1. We solve the problem numerically by using a finite element method with piecewise linear finite elements combined with either the L 2 -projection or the Carstensen interpolant of the data on the boundary. The finite element meshes for the calculations are generated by using a newest vertex bisection algorithm as described in [START_REF] Pfefferer | Numerical analysis for elliptic Neumann boundary control problems on polygonal domains[END_REF]. The discretization errors for different mesh sizes and the experimental orders of convergence are given in Tables 1-3 below for interior angles ω ∈ {3π/4, 3π/2, 355π/180}, where the discrete solutions based on the L 2 -projection and the Carstensen interpolant are denoted by y h,2 or y h,C , respectively. We note that the errors are calculated by an adaptive quadrature formula. Apparently, the results are very much in congruence with the predicted orders. The different numerical solutions y h,2 and y h,C for ω = 3π/2 and h = 1/8 are displayed in Figure 2a and Figure 2b, respectively. We see that the infinite boundary value in the origin is replaced by a finite one. If using the L 2 -projection of the data, we see also that the zero boundary values at the edge with θ = 0 are replaced by an oscillating function which is typical for this kind of regularization. By using the Carstensen interpolant as regularization of the data, this can be avoided according to the local definition of this interpolation operator. But we note that these oscillations are a feature of the former regularization approach and do not disturb the approximation order.

Extensions

Remark 5.1. The paper is written for two-dimensional domains. However most results also hold in the three-dimensional case or can simply be extended to this one. The most crucial issue is the regularity which we have used for the corresponding adjoint problem since in three dimensional domains not only corner singularities but also edge singularities need to be taken into account.

Remark 5.2. We assume f ∈ H -1 (Ω) for the discretization error estimates of Corollary 3.3. This is only for simplicity. For defining the very weak solution and the numerical solution we only need f ∈ (H 1 0 (Ω) ∩ H 1 ∆ (Ω)) ∩ Y 0h such that we could admit right hand sides from L 1 (Ω) and Dirac measures as well. However, in this case the discretization error analysis demands an adapted proof, which exceeds the scope of this paper.

Remark 5.3. We assume u ∈ L 2 (Γ) for simplicity. The case u ∈ H t (Γ), t > 0, is also of interest in the analysis of Dirichlet control problems, see [START_REF] Apel | On the regularity of the solutions of Dirichlet optimal control problems in polygonal domains[END_REF]. Proof. The interpolation error estimate (A.2) is in principle contained in [START_REF] Carstensen | Quasi-interpolation and a posteriori error analysis in finite element methods[END_REF], however, it is there an estimate on a domain such that we sketch the proof here.

For s = 0 this estimate follows from the stability property

π x (u) L 2 (ωx) ≤ c u L 2 (ωx) , (A.3)
where ω x is the support of λ x on Γ.

For s = 1 we use that π x (w) = w for all constants w such that

u -π x (u) L 2 (ωx) = (u -w) -π x (u -w) L 2 (ωx) ≤ c u -w L 2 (ωx) ≤ ch u H 1 (ωx)
via the Deny-Lions lemma, see also [START_REF] De Los Reyes | Finite element error analysis for stateconstrained optimal control of the Stokes equations[END_REF]Lemma 4.3] where the piecewise affine and Lipschitz continuous transformation of ω x to some reference domain is discussed in detail.

For a boundary edge with end points x 1 and x 2 we have

u -C h u L 2 (e) = (u -π x 1 (u))λ x 1 + (u -π x 2 (u))λ x 2 L 2 (e) ≤ c 2 i=1 u -π x i (u) L 2 (e) .
From these two estimates we obtain (A.2) in the case s = 1. In the remaining case s ∈ (0, 1) the error estimate (A.2) follows by interpolation of Sobolev spaces.

Lemma A.2. The piecewise linear Carstensen interpolant satisfies the error estimate

u -C h u H -s (Γ) ≤ ch s+t u H t (Γ) , s ∈ [0, 1], t ∈ [0, 1],
as well as the stability estimate

C h u L 2 (Γ) ≤ c u L 2 (Γ) .
Proof. The second estimate follows directly from the fact that 0 ≤ λ x ≤ 1 and the stability property (A.3). The first estimate is in principle contained in [START_REF] De Los Reyes | Finite element error analysis for stateconstrained optimal control of the Stokes equations[END_REF], however, it is there an estimate on a domain such that we sketch the proof here. First we notice that the definition of π x (u) is equivalent to (u -π x (u), λ x ) Γ = 0 and hence we have (u -π x (u), π x (ϕ)λ x ) Γ = 0 ∀ϕ ∈ H s (Γ).

With this identity and with x λ x = 1 we get

(u -C h u, ϕ) Γ = (u x λ x - x π x (u)λ x , ϕ) Γ = x (u -π x (u), ϕλ x ) Γ = x (u -π x (u), (ϕ -π x (ϕ))λ x ) Γ ≤ x u -π x (u) L 2 (ωx) ϕ -π x (ϕ) L 2 (ωx)

Figure 1 :

 1 Figure 1: Visualization of the exact solution

Figure 2 :

 2 Figure 2: Visualization of the approximate solutions with h = 1/8

  19). Lemma 2.9. Berggren's formulation (2.19), (2.20), (2.21) is equivalent to the formulation (2.6). Proof. We first assume that y ∈ L 2 (Ω) satisfies (2.6) and show (2.19)-(2.21). For any ϕ ∈ L 2 (Ω) let v be the variational solution of -∆v = ϕ defined by (2.20), hence v ∈ V and ζ := ∂ n v ∈ V Γ . Based on the formula (2.8), we obtain (2.21). With (2.6) we finally get also (2.19). Let now y satisfy (2.19)-(2.21). Since ϕ ∈ L 2 (Ω) we get from (2.20) that v ∈ V . Moreover, we obtain ζ = ∂ n v ∈ V Γ from (2.21). Hence equation (2.19) becomes

  see Lemmas 2.12 and 2.13. The finite element solution y h is now searched in Y * h := Y h * ∩ Y h and is defined in the classical way,

Table 1 :

 1 Discretization errors and experimental orders of convergence (eoc) for ω = 3π/4.

	The results can be

Table 2 :

 2 Discretization errors and experimental orders of convergence (eoc) for ω = 3π/2.

	Expected convergence rate: 1/6
	mesh size h # unknowns	y -y h,2 L 2 (Ωω)	eoc	y -y h,C L 2 (Ωω)	eoc
	0.50000	46	1.1049		1.1141
	0.25000	159	1.0693	0.04721	1.0732	0.05406
	0.12500	589	1.0491	0.02749	1.0513	0.02967
	0.06250	2265	1.0367	0.01715	1.0384	0.01782
	0.03125	8881	1.0281	0.01207	1.0296	0.01226
	0.01562	35169	1.0213	0.00956	1.0228	0.00962
	0.00781	139969	1.0154	0.00832	1.0169	0.00834
	0.00390	558465	1.0100	0.00771	1.0114	0.00772
	expected			0.00704		0.00704

Table 3 :

 3 Discretization errors and experimental orders of convergence (eoc) for ω = 355π/180, Expected convergence rate: 180/355 -1/2 = 1/142 ≈ 0.007 Lemma A.1. The piecewise linear Carstensen interpolant satisfies the error estimateϕ -C h ϕ L 2 (Γ) ≤ ch s ϕ H s (Γ) , (A.2)for ϕ ∈ H s (Γ), s ∈ [0, 1].

The work was partially supported by Deutsche Forschungsgemeinschaft, priority program 1253 and IGDK 1754.

Lemma 3.2. The finite element error estimate

holds for s = 1 2 in the convex case and for s < λ -1 2 in the non-convex case. Before we prove this lemma we can immediately imply the following final error estimate for this approach. By the triangle inequality we have

The first term is already treated in (3.4). The second term is treated in Lemma 3.2.

Corollary 3.3. The discretization error estimate

holds for s = 1 2 in the convex case and s ∈ [0, λ -1 2 ) in the non-convex case. Note that the order can be improved if the boundary datum u is more regular, see Remark 5.3.

Proof of Lemma 3.2. Let B h u h ∈ Y * h be the discrete harmonic extension defined by

The first estimate can be cited from [17, Lemma 3.2], the second follows from an inverse inequality, the third from Lemma 2.13. Now we notice that

where

Hence y 0h is the Galerkin approximation of y h 0 and therefore (since

)), we see that

improved in this more regular case: With s from Corollary 3.3 we have

The first step is again the application of Lemma 2.12, while the second can be proved in analogy to Lemma 2.13. The necessary prerequisites are already provided in Lemma A.2 and Remark A.3. Furthermore, when we check the proof of Lemma 3.2 we find that we obtain

Hence it remains to prove

in order to conclude

The estimate (5.1) is known for domains and can be proved for t ∈ [0, 1] also for the boundary Γ by using the inverse inequality and the approximation properties of u h and the Ritz projection. To this end recall that Γ j , j = 1, . . . , N , are the boundary segments of Γ, and let P j h be the Ritz projection on Y ∂ h | Γ j . Then we have

i. e., we have (5.1) for t = 1. Since we proved (5.1) for t = 0 in Lemma 2.13 we get the desired result by interpolation in Sobolev spaces.

Appendix: Error estimates for the Carstensen interpolant

Recall from Subsection 2.5 that the piecewise linear Carstensen interpolant is defined via

where λ x is the standard hat function related to x.

where ω x is again the support of λ x on Γ. With similar arguments as in the proof of Lemma A.1 we conclude

By the definition of the negative norm,

we obtain the assertion of the lemma.

Remark A.3. Note that this error estimate holds also for the L 2 -projection,

It can be proved similarly by using (u -Π h u, ϕ) Γ = (u -Π h u, ϕ -Π h ϕ) Γ .

Acknowledgement The authors thank Markus Melenk and Christian Simader for helpful discussions.


Susanne Fuchs
email: fuchs@zas.gwz-berlin.de

Pascal Perrier

Christian Geng

Christine Mooshammer

What role does the palate play in speech motor control? Insights from tongue kinematics for German alveolar obstruents

By means of simultaneous EMMA and EPG recordings, we investigated tongue tip kinematics and tongue palate contact patterns for four German speakers in order to compare production strategies of alveolar stops with fricatives. For alveolar stops versus fricatives, two different control strategies were hypothesized: a target above the contact location for alveolar stops resulting in a collision of the tongue tip at the palate as opposed to a precise positioning of the tongue at the lateral margins at the palate for alveolar fricatives.

Additionally, we suspected differences between stops and fricatives with respect to anterior and lateral palate contacts and their influence on tongue kinematics. Results of this study confirmed two different control strategies for alveolar stops and fricatives by means of significant differences in movement amplitude, velocity, and duration of the closing gesture, the amplitude of deceleration peaks, tongue tip movement during acoustically defined closure and constriction, and maximal anterior contact during closure.

Additionally, results for speaker dependent mechanisms were related to the subject's coronal palatal shape.

Introduction

A particular characteristic of speech production compared to other human motor systems is that the tongue moves in a narrow space delimited by the pharyngeal walls, cheeks, teeth, lips and palate. Hence, even if the basic principles that govern arm or limb motor control might also apply to speech, it seems reasonable to assume that in speech production specific strategies are used to deal with these additional constraints. The focus of this paper is on the upper extent the tongue's possible movement, the palate. Since tongue movements for oral consonants are most of the time produced in the presence of palatal contacts, it is hypothesized that: 1) the palate passively influences kinematic properties of speech movements, and 2) the tongue's action at the palate is taken into account in the speech motor control process in terms of limiting the degrees of freedom for tongue movement and supporting the tongue's shaping. As far as consonant production is concerned, [START_REF] Stone | Toward a model of three-dimensional tongue movements[END_REF] suggested that some tongue shapes, particularly those in the production of alveolar fricatives, could not be produced by a free-standing tongue position, i.e., with the palate as a spatial reference for the tongue; Stone also suggests that the palate can be used to provide tactile feedback information in order to learn and/or control specific tongue shapes. Additional evidence for the role of the palate was provided by perturbation studies. [START_REF] Hamlet | Compensatory alveolar consonant production induced by wearing a dental prosthesis[END_REF] Stone recorded tongue-palate contact patterns for ten speakers wearing a dental prosthesis of 4 mm thickness in the alveolar region. Data were observed in two sessions: immediately after the insertion of the prosthesis, when subjects were still unfamiliar with it, and two weeks after adaptation. Results from the first session provided evidence for tongue overshoot in /s, z, t, d, n/. Results for the second session showed compensatory movements with a similar variability as in the first session. The authors suggested that:

The way sensory feedback is used in learning a compensatory form of speech is qualitatively similar to its use routinely. That is, all sensory information (from audition, proprioception, touch etc.) by which the status of the oral environment and relative positions of the articulators can be assessed, is continuously being monitored, and used in planning the speech act (Hamlet and Stone 1978, p.246). [START_REF] Baum | The development of speech adaptation to an artificial palate[END_REF] investigated acoustic and perceptual changes for /s/ production after structural modifications of the alveolar ridge using an artificial palate. Seven subjects were analyzed. They found that /s/ was acoustically and perceptually highly susceptible to perturbing effects in the beginning of the session. After one hour's practice, a significant improvement of acoustic and perceptual characteristics was observed, although compensation was not complete. The authors propose that an adaptation to articulatory perturbations could be achieved in a relatively short time period due to a recalibration of the speech motor control system taking into account the new spatial reference at the alveolar ridge. [START_REF] Honda | Compensatory responses of articulators to unexpected perturbations of the palate shape[END_REF] and [START_REF] Honda | Effects of tactile and auditory feedback on compensatory articulatory response to an unexpected palatal perturbation[END_REF] used another experimental paradigm, dynamic (time varying) perturbations. They examined compensatory tongue movement for unexpected perturbations due to an inflatable artificial palate whose thickness could suddenly be increased or reduced; this technique was combined with or without auditory masking. Their findings suggest that incomplete compensatory movements are actively and immediately produced due to tactile feedback information, and that auditory feedback is used to complete precise articulatory adjustment, but with a longer time delay in comparison to tactile feedback.

The influence of the tongue's action at the palate on tongue kinematics was discussed in [START_REF] Hoole | Theoretische und methodische Grundlagen der Artikulationsanalyse in der experimentellen Phonetik[END_REF] who found for German /aCa/-sequences differences with respect to the deceleration peak of the closing gesture in the following order: /t/>/d/>/n/>/s/. He assumed that in alveolar stop production the relatively high deceleration peak of the tongue tip is mainly passive and a result of tongue's action at the palate whereas in alveolar fricatives the deceleration phase is actively controlled. [START_REF] Fuchs | The role of the palate in tongue kinematics: an experimental assessment in VC sequences from EPG and EMMA data[END_REF] showed that the onset of oral closure, defined by EPG data, coincided with the time point of the high deceleration peak for the tongue tip sensor. It was proposed that a collision of the tongue tip against the palate could be responsible for the abrupt deceleration and the strength of the deceleration peak.

These findings support the hypothesis of [START_REF] Löfqvist | Lip and jaw kinematics in bilabial stop consonant production[END_REF] for bilabial stops that articulatory gestures could be directed toward a target that is beyond the contact location.

Thus, coronal stops are not controlled via a fine positioning of the tongue at the palate, but by the use of the palate as a reference which automatically blocks the tongue's movement at the required location in the vocal tract. In terms of stability and simplicity, such a control strategy seems to be extremely efficient in comparison with the control of a fine positioning of the tongue at the palate. The effects of target variations in stop production on tongue kinematics (a precise positioning of the tongue without the palate and a collision between tongue and palate) were modeled in [START_REF] Perrier | Influences of tongue biomechanics on velar stop consonants: A modeling study[END_REF] using a 2D biomechanical tongue model. They showed that the tongue's action at the palate combined with the relative articulatory target positions of the surrounding phonemes could explain some of the direction of the anticlockwise trajectories (so called loops) in /VkV/-sequences.

In contrast to stops, [START_REF] Mooshammer | Coordination of lingual and mandibular gestures for different manners of articulation[END_REF] observed that fricatives have a more precise positioning of tongue and jaw. Fricatives are also often described as being difficult to produce. They are acquired later in the children's speech development in comparison to less complex sounds, cause most problems in patients who suffer from sensori-motor coordination impairments and they seem to occur more frequently in tongue slips, especially regarding place of articulation [START_REF] Boomer | Slips of the tongue[END_REF].

Previous studies have mainly focused on tongue kinematics by means of EMMA data and have inferred the corresponding tongue palate contact patterns. We will take the advantage of simultaneous EMMA and EPG recordings. The aims of our study are twofold. First we will compare alveolar stops with fricatives, because fricatives show both a fine positioning of the tongue at the lateral margins of the palate as well as tongue grooving which is not required in stop production. Thus, we suppose two different control strategies for stops versus fricatives.

Second, the tongue's action at the anterior and lateral parts of the palate were investigated separately in order to interpret their potential influence on tongue kinematics. Lateral contacts could potentially cause damping, and they might explain speaker dependent differences in tongue kinematics based on the shape of the palate (Tiede 1998, [START_REF] Mooshammer | An EMMA and EPG study on token-to-token variability[END_REF]).

Method

Tongue tip movements and tongue-palate contact patterns from four German speakers (all colleagues from the lab), three male (CG, DF, JD) and one female (SF, the first author of this paper), were recorded simultaneously by means of EPG (Reading EPG3) and EMMA (AG100, Carstens Medizinelektronik) systems. Tongue tip movement was associated with the movement of a sensor placed mid-sagittally approximately 1 cm behind the tip. Two sensors served as reference points to compensate for helmet movements, one at the bridge of the nose and one at the upper incisors. All other sensors at the tongue and the jaw were not taken into account here. Speech signals were synchronously recorded on DAT. The sampling frequencies were 16 kHz for the acoustic data, 100 Hz for EPG and 200 Hz for EMMA data respectively. The speech material consisted of the nonsense words ge/C1VC2/e and ge/C1VC2/, where C was always either /t/ or /z/. The phonologically voiced /z/ has been used instead of /s/ since it coincides most frequently with a preceding tense vowel in the German inventory whereas /s/ occurs with lax vowels. The focus of this paper is on C2 of each target word, following one of the stressed vowels // or /u/. Different positions of the consonant were considered: word medial C (henceforth med) in e.g., [f?!sç@9s?] corresponding to orthographic 'getahte', and word final position C (henceforth fin) in e.g., [f?!sç@9s], orthographically 'getaht'. All target words were embedded in the carrier phrase "Ich habe geC1VC2e nicht geC1VC2 erwähnt." (= I have geC1VC2e not geC1VC2 said, i.e. I said geC1VC2e not geC1VC2.) with pitch-accents on C1VC2 in both cases to convey contrastive sentence stress. A neutralization process of word final /z/ to the voiceless allophone takes place in German, i.e. word final /z/ was in most cases realized as [s]. As will be seen in the results section, supralaryngeal production mechanisms of word medial phonological /z/ (produced as [z]) and word final /z/ (produced as [s]) did not differ most of the time so that voicing distinction was not likely to be relevant for the interpretation of the tongue kinematics in our data. Thus, we maintained the phoneme description /z/.

Each sentence was repeated ten times in a randomized order. The tongue tip closing gesture from the stressed V to C was measured. The on-and offsets of the closing gesture were defined using a 20 per cent threshold criterion (see figure 1), i.e. closing gesture onset corresponds to the point where the tangential velocity signal of the tongue tip reached 20 per cent by traveling from the velocity minimum to the peak. Closing gesture offset corresponds to the time landmark where the velocity signal reached 80 per cent by traveling from the velocity peak to the next velocity minimum. The advantages of the 20 per cent threshold criterion are described in [START_REF] Kroos | Phonetic evidence for the phonological status of the tense-lax distinction[END_REF].

Insert Figure 1 around here Starting from the landmarks defining the closing gesture, the velocity peak, the duration of the closing gesture and the movement amplitude were calculated. The latter was calculated as the integral of the tangential velocity signal from onset to offset of the closing gesture which corresponds to the area beyond the velocity signal. The measurement is generally more precise in comparison to calculations of the Euclidian distance, since it takes into account the whole curvature of the trajectory.

In order to compare stop production with fricative production in relation to previous findings in [START_REF] Fuchs | The role of the palate in tongue kinematics: an experimental assessment in VC sequences from EPG and EMMA data[END_REF] and [START_REF] Hoole | Theoretische und methodische Grundlagen der Artikulationsanalyse in der experimentellen Phonetik[END_REF], the acceleration and deceleration peaks of the closing gesture (VC-sequence) as well as of the preceding opening gesture (CV-sequence) were labeled. Peaks in the preceding opening gesture were taken into account in order to have some references for the magnitude of the deceleration peaks.

Two acoustic landmarks were used for analyzing the obstruents. For the stops, the closure and burst were labeled and the closure onset was marked at the offset of F2 in the preceding vowel. The acoustic boundaries for the fricatives were marked at the on-and offset of high frequency noise.

From each EPG pattern, the percentage of contact in the anterior region (4 most front rows) = ANT, in the posterior region (4 most back rows) = POST and in the lateral region (2 most peripheral columns on the left and right side of the palate) = LAT were computed in order to quantify tongue palate contact in time. Onsets of lateral and anterior contacts of the closing gesture were labeled as the landmarks where the relevant contact increased after the preceding vowel. The increase of contact was most of the time rather abrupt, particularly for stops which showed nearly no contact during the end of the vowel and a sudden rise in the following pattern for closure onset. Additionally, the landmarks of maximal lateral and anterior contact during oral closure/constriction were marked. In order to get an impression of subject's palatal shape, for each of the 62 electrodes placed in the individual artificial palates, the x, y, and z coordinates were measured using a caliper (for a detailed description of the method see [START_REF] Fitzpatrick | Estimating lingual constriction location in high vowels: a comparison of EMA-and EPG-based measures[END_REF].

Insert Figure 2

Results

Closing gesture amplitude and its duration

A number of analyses of variance were calculated using the General Linear Model in SPSS 11.5. Closing gesture amplitude, its duration and velocity peak were the dependent variables and consonant (/z/ versus /t/) was the independent variable. The data were split by vowel (/a/ versus /u/) and position in word (med versus fin). The results showed that movement amplitudes in /t/ were significantly larger (see figure 3) for all subjects in all conditions (all p < 0.001, except JD med /a/-context p < 0.05 at a 95 per cent confidence interval). The same was true for the velocity peaks (all highly significant).

Insert Figure 3 around here Similar findings have been presented by [START_REF] Geumann | Invariance and variability in articulation and acoustics of natural perturbed speech[END_REF] who showed a smaller movement amplitude due to a lower tongue tip target position in /s/ and a higher tongue target position in /t/. The smaller movement amplitude for the fricative may also be a result of tongue grooving as described e.g. in [START_REF] Narayanan | An articulatory study of fricative consonants using magnetic resonance imaging[END_REF], and because the tongue tip sensor was placed within the mid-sagittal groove. Tongue grooving may also be confirmed by means of maximal anterior contact during oral closure for the stops in comparison to the contacts during oral constriction for the fricatives. Results for stop production showed generally significant more per cent of maximal anterior contact than in fricative production (all /a/contexts p < 0.001, except DF fin p < 0.05, /u/-context med CG, SF p < 0.001, JD p < 0.05, fin CG, SF p < 0.001, JD p < 0.05). The smaller amount of anterior contact in /z/ should be due to the missing palatal contacts in the tongue groove.

If alveolar fricatives and stops were controlled in the same way, we would expect movement duration to be similar for both types of consonants. We would also expect larger movement amplitudes to coincide with greater velocity peaks and smaller movement amplitudes with smaller velocity peaks. Such results were generally found for the alveolars realized in /u/context (see figure 4). However, in all speakers' production of the /a/-context, the movement amplitude and the velocity peak were larger in /t/ compared to /z/, but the closing gesture duration was shorter for /t/ than /z/ (med: CG p < 0.001, DF/JD/SF p < 0.01; fin: CG/DF p < 0.001, SF p < 0.01, JD p < 0.05).

The shorter duration observed for /t/ in the /a/-context gives a quantitative information about the consequence of the differences observed in peak velocity and movement amplitude between fricatives and stops. If the movement amplitude is large enough, the tongue can move in a shorter time interval for stops as compared to fricatives. In the case of the /u/context these differences cannot be observed since the movement amplitude was too small.

Insert Figure 4 around here

Deceleration of tongue tip movement and movement during closure/constriction

Two separate measurements were made in order to further assess strategies used in stop and fricative production: first, acceleration and deceleration peaks in the opening and closing gestures were compared and second, tongue tip movements (movement amplitude, vertical and horizontal movements) during the oral closure or constriction were analyzed.

As predicted from the hypothesis that the tongue's collision at the palate strongly affects movement kinematics for stops, a consistently higher deceleration peak during the closing gesture was found in comparison to the opening gesture for /t/ (see figure 5). This was also true for /z/ (except for SF). Additionally, the relative differences between deceleration peaks in the opening and closing gestures were always greater for /t/ than /z/, e.g. for CG -the deceleration peak in the closing gesture for /t/ was 2.2 times greater than in the opening gesture whereas in /z/ it was only 1.4 times greater. Acceleration peaks did not exhibit comparable relations and varied considerably between speakers. In order to rule out a possible influence of the larger movement amplitude in /t/ on the deceleration peak (this could come due to the fact that larger amplitudes are associated with larger peak velocities, and thus, they are necessarily associated with deceleration, because deceleration is the velocity's derivative) a covariance analysis was calculated with deceleration peak in the closing gesture as the dependent variable, movement amplitude as the covariate, and consonant (/t/ versus /z/) as the independent variable. Data were split by position (med versus fin) and vowel context (/a/ versus /u/). Significant differences for the deceleration peaks between /t/ and /z/ were found in all /a/-context (med: CG, JD p < 0.001, DF p < 0.01, SF p < 0.05, fin: CG, DF, SF p < 0.001, JD p < 0.001). In /u/-context differences were less consistent (med: JD, SF p < 0.001, DF p < 0.01, fin: SF p < 0.001, DF p < 0.05). These findings suggest that the tongue's movement does not stop due to an active deceleration in alveolar stop production, but due to an impact at the palate. It is in general agreement with Hoole's proposal (1996) and Löfqvist's and Gracco's suggestion (1997) for a planned target position beyond the contact location. patterns observed for this speaker). Forward movements during /t/ closure (i.e. the tongue tip slides along the palate) were more pronounced in /u/ than in /a/-context, confirming the idea that a redirection of the movement depends on the planning of the target and on the angle of incidence. [START_REF] Brunner | Correlation between angle of incidence and sliding patterns of the tongue along the palate in Korean velar stops[END_REF] suggested that the greater the size of the angle of incidence, determined by the target positions of the surrounding environment, the larger the movement amplitude of the sliding tongue. Although the exact angle of incidence was not calculated here, it seems reasonable to assume that the angle of incidence in /u/-context is more obtuse than in /a/-context since the tongue moves from a back high vowel towards the front.

Duration from onset of lateral contacts to velocity peak

In a previous pilot study of one subject, we observed that the velocity peak in /t/ was well synchronized with an increase in lateral contacts, corresponding to the onset of stop occlusion. This observation raised the question whether for stop production the deceleration phase could be a result of tongue's action at the palate only. In order to begin to test the reliability of these preliminary findings, the temporal difference between the velocity peak and the onset of lateral contact was calculated. In the /u/-context (all speakers, all conditions) lateral contacts started already during the vowel, i.e. they preceded the velocity peak. Tongue and palate were always in contact with each other since /u/ is a high vowel. Results obtained for the /a/-context are shown in figure 8. They show considerably inter-speaker variation.

Insert Figure 8 around here Differences between /z/ and /t/ are relatively small keeping in mind that the frame rate of EPG is 100 Hz. However, results clearly reveal specific patterns (/a/-context) for speakers with a dome shaped palate versus speakers with a flat palatal shape. For the latter (JD and SF), tongue contact with the lateral margins of the palate occurred already during the preceding /a/ (between 15-30 per cent LAT) and increases during the consonant. This type of pattern did not occur in speakers with a dome shaped palate (CG and DF). Their results showed roughly a simultaneous timing of the velocity peak with the onset of lateral contacts.

In contrast, for speakers with a flat palatal shape (JD and SF), the tongue is already in contact with the palate during the vowel. One consequence of this seems to be a weaker impact of the tongue at the palate during stop production, as suggested by comparing the magnitude of the deceleration peaks during /t/ for the two groups of speaker. Those with a flat palate had lower deceleration peaks, suggesting that contact with the palate caused less damping of the tongue movement than it is the case for the speakers with more dome palate shapes. In this way, the variations in the tongue palate contact patterns during the preceding vowel affected tongue kinematics in different ways. Since lateral contacts were already produced during the acceleration phase, they cannot induce the deceleration phase alone. However, a combination of lateral and anterior contacts may cause the beginning of the deceleration, at least for stop production.

Contrary to our expectation that lateral contacts start earlier with respect to anterior contacts in /z/ production, we did not find reliable differences between /t/ and /z/. If differences occurred, then they were often 10 ms and therefore close to the reliability range of EPG. This means that either lateral contacts always co-occured with anterior contacts in alveolar obstruent production or that differences between the onsets of such contacts were very tiny and should be investigated by means of a higher time resolution technique.

Discussion

We investigated tongue tip kinematics and tongue palate contact patterns by means of simultaneous EMMA and EPG recordings for four German speakers in order to compare production strategies of alveolar stops with fricatives. We assumed that the central nervous system incorporates the palate in speech motor control and that tongue palate contact patterns influence the kinematics of the tongue movement. In stop production, the alveolar target is planned beyond the contact location (see [START_REF] Löfqvist | Lip and jaw kinematics in bilabial stop consonant production[END_REF]. In terms of stability and simplicity, such a control strategy would be efficient, since it can tolerate quite a large amount of variability in the motor commands without endangering the final tongue positioning and the corresponding auditory goal. The tongue's collision at the palate guarantees an air tight seal corresponding to a silent closure phase in the acoustics. By contrast, fricatives are planned with a target at the lateral margins of the palate in order to create a narrow mid-sagittal constriction as a prerequisite for the production of turbulent noise, the most prominent acoustic characteristic of fricatives. To do so, a more precise tongue and jaw position is required.

Results for both strategies were evident in tongue tip kinematics and tongue palate contact patterns. The larger deceleration peak in /t/ production provides evidence for a collision of the tongue tip with the palate (in agreement with [START_REF] Hoole | Theoretische und methodische Grundlagen der Artikulationsanalyse in der experimentellen Phonetik[END_REF][START_REF] Fuchs | The role of the palate in tongue kinematics: an experimental assessment in VC sequences from EPG and EMMA data[END_REF].

Additionally, the movement amplitude and the velocity peak of the closing gesture were larger in /t/ whereas the duration of the closing gesture was significantly shorter, at least in the /a/-context. Movement amplitudes for /z/ were smaller, slower, and were produced with a longer duration. There was a significantly greater percentage of maximal anterior contact during closure or constriction for /t/ than for /z/, supporting two ideas: First, alveolar fricatives are produced with tongue grooving and have therefore less contact with the palate; second, an impact always produces a high amount of contact whereas fine positioning can potentially produce either a small or a large amount of contact.

By observing the movement amplitude during the acoustically defined oral closure for stops and the oral constriction for fricatives, a larger tongue tip movement amplitude was found for /t/ than for /z/. The tip generally moves in a forward direction which is particularly pronounced in the /u/-context and confirms earlier work from [START_REF] Brunner | Correlation between angle of incidence and sliding patterns of the tongue along the palate in Korean velar stops[END_REF]. In /z/ production, less movement was produced during oral constriction, probably because of the more precise tongue positioning that is required in fricatives than in stops. We suppose that the different control strategies are used in order to create the appropriate auditory goals, i.e. a silent closure interval for stops and turbulent noise for fricatives.

Concerning the role of lateral tongue palate contacts, differences were less consistent for the two consonants. However, speaker dependent variations were found and interpreted with respect to differences in palate shape. Speakers with a flat palate shape exhibit lateral contact already during /a/-production whereas for speakers with a dome shaped palate such contact did not occur. As far as the /a/-context is concerned, we suggest that for speakers with a flat palate shape the impact between tongue and palate was less pronounced since the tongue was already in contact with the palate during the preceding vowel. A previous study of the tokento-token variability of all German vowels [START_REF] Mooshammer | An EMMA and EPG study on token-to-token variability[END_REF]), including the same male speakers recorded here, provided evidence for a smaller variability of tongue positioning in general, comparing JD with CG and DF. We interpreted this result partly due to differences in palate shape, since small tongue movements have larger effects on the area function for speakers with a flat palate. Thus, speakers with a flat palate may have learnt a more precise tongue positioning. Whether this also holds true for consonant production and the extent to which differences in the precision of vowel production affect consonant production need further investigation. [START_REF] Johnson | Individual differences in vowel production[END_REF] proposed that palatal doming affects articulatory organization, in particular tongue and jaw co-ordination. In a subsequent investigation, we will also observe jaw movements, since preliminary inspections showed differences in vertical jaw movements for the two groups of speakers with smaller movement amplitudes for speakers with the flat palate (JD, SF) and larger amplitudes for the speaker with a dome shaped palate (CG, DF).

List of figures

Figure 1: Articulatory and acoustic labeling criteria from CG's medial /t/ production in /a/context: 1 st track: audio signal with f2 offset and the burst, 2 nd track: tongue tip vertical movement in cm, 3 rd track: tongue tip tangential velocity signal in cm/s with closing gesture onset and offset, 4 th track: tongue tip tangential acceleration signal in cm/s 2 with acceleration (acc) and deceleration peaks (dec)

  Figure 2 exhibits the palatal shapes for all subjects based on the coordinates of the 62

Figure 2 :

 2 Figure 2: Coronal shapes of the artificial palates (anterior to posterior view for four speakers), positions of the 62 electrodes are marked by the open circles, lines connect open circles in the relevant row and column, x-axis (mm): left-right axis, y-axis (mm): palatal height

  Figure 2

Acknowledgements

This work was supported by a grant from the German Research Council (DFG) GWZ 4/8-1, P.1. We thank Jonathan Harrington and two anonymous reviewers for their helpful comments on an earlier version of this manuscript, our subjects, and Jörg Dreyer for technical support during the recordings. Tiede, M. (1998). Palate shape effects on characteristic vowel tongue postures. Journal of the Acoustical Society of America, 103, 2893.

Additional evidence was derived from the results concerning tongue tip movement during acoustically defined closure and constriction. If a target is planned beyond a contact location, the tongue's movement is reoriented due to the collision at the palate. The reorientation of the movement depends on the angle of incidence the tongue collides at the palate. For example, if the tongue moves against the palate in an obtuse angle coming from the back it should be redirected towards the front and sliding along the palate. If the angle is around 90 degrees the tongue should have the greatest impact on the palate since the force cannot be redirected (for calculation of the angle of incidence by means of tongue kinematics in looping patterns see [START_REF] Brunner | Correlation between angle of incidence and sliding patterns of the tongue along the palate in Korean velar stops[END_REF]). For a target planned at the contact location, less damping effects should occur due to the missing impact, and less movement should be produced in order to place the articulator precisely at the target location.

Insert Figure 6 around here Figure 6 shows that there were larger movement amplitudes during /t/ closure compared to /z/ constriction in all cases for the /u/-context (highly significant). In the /a/-context, significant differences between /t/ and /z/ occurred particularly in the word medial position (DF, JD, SF).

In word final position a /t/-/z/ distinction was only found for DF. Results from subject CG showed generally the reverse pattern in /a/-context with larger amplitudes for /z/ than for /t/.

Insert Figure 7 around here The larger movement amplitudes during oral closure for /t/ were mainly a result of a forward movement of the tongue tip as shown in figure 7 (except for CG in /a/ context for whom differences in amplitude corresponded to vertical tongue movement, which explains the