全球半導(dǎo)體晶體生長建模著名商業(yè)軟件femag--塊晶體生長工業(yè)級數(shù)值模擬

1、,全球半導(dǎo)體晶體生長建模著名商業(yè)軟件FEMAGNumerical Simulationof Bulk Crystal Growthfor Industrial Application,François Dupret1,2, Roman Rolinsky2, Brieuc Delsaute2, Rajesh Ramaya2, Nathalie Van den Bogaert2 1 Université ca

2、tholique de Louvain, Louvain-la-Neuve, Belgium 2 FEMAGSoft S.A. company, Louvain-la-Neuve, Belgium,,FEMAGSoft © 2013,How to improve the growth process in terms of: crystal quality ? process yield ?

3、 energy consumption ? production rate ?,Introduction,Main difficulties:,,FEMAGSoft © 2013,Multi-physics: heat and mass transport in the melt and the gas, turbulence, radiation transfer, etc., all interact and stron

4、gly affect species incorporation and defect formation in the crystal,Multiple space scales: sharp diffusive, viscous, radiative and thermal boundary layers are present in the melt and the gas, together with complex defec

5、t boundary layers in the crystal,Multiple time scales: typically the growth process is very slow while the melt flow is governed by much shorter time constants,Introduction (cont’d),FEMAGSoft © 2013,Solving these pr

6、oblems requires …,To resort to appropriate and up-to-date numerical simulation techniques to couple and solve these models→ quasi-steady and dynamic models,To develop a sound physical model for each separate effect→ gl

7、obal and time-dependent modeling of heat transfer, turbulence modeling, defect modeling, …,,Introduction (cont’d),,FEMAGSoft © 2013,Principal objective has been to complete the platform,FEMAG-2 → FEMAG-3 software ge

8、neration transition taking place from 2008-2009,General objective of FEMAGSoft,→ strongly improved platform in terms of computation time, memory, etc.,Introduction (cont’d),FEMAGSoft © 2013,b) Time-dependent modelin

9、g: use of various simulation modes(ex: quasi-steady, quasi-dynamic, inverse or direct dynamic models in Cz growth),c) FEM discretization: use of 2D, Spectral 3D, Cartesian 3D,… models(high geometrical flexibility, si

10、mple assembling technique),a) Global modeling: subdivision of the furnace into “macro-elements”(solid or liquid constituents, radiation enclosures, “cement” elements...),d) Geometrical modeling: to accurately handle st

11、rongly deforming bodies and interface and well-capture all the boundary layers,e) Solution technique: coupled Newton-Raphson iterations by use of a highly effective linear solver,,Introduction (cont’d),FEMAG software dev

12、elopment strategy,1. Numerical strategy,,FEMAGSoft © 2013,Numerical strategy,FEMAGSoft © 2013,Quasi-steadythermal equilibriumadapted heater power to get the prescribed crystal diameterheat source on the soli

13、dification front in proportion to the pull rate,Inverse dynamicadapted heater power to grow the prescribed crystal shapeeffect of pull rate and solid-liquid interface deformation on the solidification heat,Direct dynam

14、iccalculated crystal shapeprecribed heater power historyeffect of pull rate and solid-liquid interface deformation on the solidification heat,Time dependent,Quasi-steady,Quasi-dynamicfrozen geometry (except the solid

15、-liquid interface)adapted heater power to get the prescribed crystal diametereffect of pull rate and solid-liquid interface deformation on the solidification heat,Different simulation modes,1. Numerical strategy (cont’

16、d),,,Global temperaturefield,Streamfunction,FEMAGSoft © 2013,Inverse QS and TD simulation of the growth of a 300 mm silicon crystal,Analysis of conical growthand shouldering stagesm = 8.225 10-4 kg/m.sWc= 3.82

17、rpm (0.4 s-1)Ws= -3.82 rpm (-0.4 s-1)Vpul = 1.8 cm/h (5. 10-6 m/s),1. Numerical strategy (cont’d),,FEMAGSoft © 2013,Temperaturefield,Streamfunction,1. Numerical strategy (cont’d),,FEMAGSoft © 2013,FEMAG-1 t

18、ime-dependent simulation of Czochralski Ge growth,Inverse dynamic simulation (imposed crystal shape, calculated heater power): power oscillations resulting from inverse modeling, and smoothed power,1. Numerical strategy

19、(cont’d),,FEMAGSoft © 2013,Direct dynamic simulation (imposed stepwise decrease of heater power, calculated crystal shape): evolution of the temperature field,FEMAG-1 time-dependent simulation of Czochralski Ge grow

20、th,1. Numerical strategy (cont’d),,Inverse dynamicoften more reliable than quasi-steady modelhighly attractive to predict crystal quality,Quasi-steadyfrequently usedcheap, but not always validdoes not allow crystal

21、quality prediction,Direct dynamicsimulation of the system response to perturbations of the input parametersvery useful for controller design,Quasi-dynamicmay capture the detailed system dynamics at various stagesvery

22、 useful for controller design,FEMAGSoft © 2013,Different simulation techniques,1. Numerical strategy (cont’d),,,Typical FEMAG-CZ global unstructured mesh,FEMAGSoft © 2013,,Heat shield,,components,1. Numerical s

23、trategy (cont’d),FEMAG-CZ global unstructured mesh deformation,FEMAGSoft © 2013,1. Numerical strategy (cont’d),,,Detail,FEMAGSoft © 2013,1. Numerical strategy (cont’d),,Use of special CAGD techniques for menisc

24、us calculation during tail-end stage : Secondary mesh (constraining loci) and deformed melt- crystal and melt-gas interfaces (right) Generated 1D and 2D meshes (left),FEMAGSoft © 2013,1. Numerical strategy (con

25、t’d),,Smooth effect of remeshing,FEMAGSoft © 2013,1. Numerical strategy (cont’d),,,,,,,,,,,,,,,,,,t0,t1,t2,t3,t4,t5,t6,time,,,t7,,,Cone growth,Body growth,Tail-end stage,Mesh re-generation,Mesh deformation,Mesh gene

26、ration,,,,,General geometrical strategy,,FEMAGSoft © 2013,1. Numerical strategy (cont’d),,FEMAGSoft © 2013,The FEMAG-3 platform,FEMAG-1,FEMAG-2,FEMAG-3,1990,2008,1984,1. Numerical strategy (cont’d),,,FEMAGSoft

27、© 2013,Well-organized architecture:,Transition from previous FEMAG-2 generation almost complete,Very fast solver (favourably competes with existing open source software),Basic features,1. Numerical strategy (cont’d)

28、,-Level 0: basic linear algebra operations (LAM),-Level 1: advanced linear algebra,-Level 2: Finite Element discretization modules,-Level 3: geometrical modules,-Level 4: non-linear solver,-Level 5: time-dependent

29、solver,- considerable saving of computer resources (computer time, memory, etc.).,,FEMAGSoft © 2013,As a result:,- considerable saving of development time.,1. Numerical strategy (cont’d),Global Simulation of Czochra

30、lski Silicon Growth under the Effectof a Transverse Magnetic Field,FEMAGSoft © 2013,,2. Cz Si growth under a TMF,,MHD boundary layers: Hartmann layers,,FEMAGSoft © 2013,The Hartmann layers develop along the su

31、rfaces where the normal component of the magnetic field is non-negligible.,An order of magnitude of the thickness of a Hartmann layer is dH = L Ha-1 (L = Rs or Rc).Typically dH = 0.05 - 0.08 mm in industrial fu

32、rnaces.,,2. Cz Si growth under a TMF (cont’d),,Transverse magnetic fields:FLET method,Method: Fourier decomposition of all fields (velocity, pressure, temperature),Hypothesis:,Principal issue: which and how many modes

33、?,Objective: global, quasi-steady or time-dependentcalculations at a reasonable cost,,,,FEMAGSoft © 2013,,2. Cz Si growth under a TMF (cont’d),Fourier Limited Expansion Technique (FLET),FEMAGSoft © 2013,The

34、field variables are expanded as Fourier series in the azimuthal (q) direction:,,2. Cz Si growth under a TMF (cont’d),Fourier Limited Expansion Technique (FLET),FEMAGSoft © 2013,and while the number of modes M is a

35、s small as possible without loss of accuracy (spectral convergence).,,while the different mode coefficients:,are 2D Finite Element functions of the meridional coordinates (r, z),,,,This results in a system whose size is

36、 that of the 2D system multiplied by the number of Fourier modes considered and hence in a dramatic system size reduction.,2. Cz Si growth under a TMF (cont’d),,FEMAGSoft © 2013,Radiation transfer,Coupling between 3

37、D and 2D axisymmetric heat transfer across a radiative enclosure,the viewed and hidden parts are calculated as axisymmetric or, equivalently, each surface of the enclosure is viewed as axisymmetric from the other surf

38、aces,Main modeling hypothesis:,,generally 3D components are rotating with respect to the 2D environment 3D components mostly view 2D components because of the presence of heat shields,This hypothesis is satisfactory

39、 because:,,2. Cz Si growth under a TMF (cont’d),FEMAGSoft © 2013,Flow and global heat transfer in a silicon Cz pullerunder the effect of a TMF(quasi-steady simulation),,2. Cz Si growth under a TMF (cont’d),FEMAGSo

40、ft © 2013,Top: top view of the velocity magnitude and streamlines.Bottom: velocity field magnitude and cross-section showing a sharp Hartmann layer along the melt-crucible interface.,Growth of a 300 mm diameter S

41、i crystal under the effect of a 0.5 T TMF.,,2. Cz Si growth under a TMF (cont’d),FEMAGSoft © 2013,In a typical TMF configuration extremely thin Hartmann layersof 50-80 mm develop,Growth of a 300 mm diameter Si

42、 crystal under the effect of a 0.5 T TMF.,,2. Cz Si growth under a TMF (cont’d),FEMAGSoft © 2013,Detail of the deforming Boundary Layer Mesh (BLM) used,Flow and global heat transfer in a silicon Cz pullerunder th

43、e effect of a TMF,,2. Cz Si growth under a TMF (cont’d),FEMAGSoft © 2013,Hartmann boundary layers along the melt-crystal and melt-crucible interfaces and associated boundary layer meshes.Strong Hartmann backflows

44、develop.,Growth of a 300 mm diameter Si crystal under the effect of a 0.5 T TMF.,,2. Cz Si growth under a TMF (cont’d),FEMAGSoft © 2013,,2. Cz Si growth under a TMF (cont’d),Czochralski growth of a silicon crysta

45、lunder a 500 mT horizontal magnetic field,Global simulation the growth of 300 mm and 400 mm crystals,FEMAGSoft © 2013,,2. Cz Si growth under a TMF (cont’d),Growth of a 300 mm crystal under a 500 mT TMF,Le

46、ft: melt surfaceRight: meridional cross-sections parallel and perpendicular to the magnetic field,Top: velocity fieldBottom: temperature field,FEMAGSoft © 2013,,2. Cz Si growth under a TMF (cont’d),Growth

47、of a 400 mm crystal under a 500 mT TMF,Left: melt surfaceRight: meridional cross-sections parallel and perpendicular to the magnetic field,Top: velocity fieldBottom: temperature field,FEMAGSoft © 2013,

48、,Czochralski growth of a silicon crystalunder a 3000 or 5000 G horizontal magnetic field,Global simulation of the growth of 300 mm and 400 mm crystals,2. Cz Si growth under a TMF (cont’d),FEMAGSoft © 2013,,Ge

49、ometry and operating conditions,2. Cz Si growth under a TMF (cont’d),FEMAGSoft © 2013,,Description of the simulations,2. Cz Si growth under a TMF (cont’d),FEMAGSoft © 2013,,Crystal diameter: 400 mmPulling rate

50、: 0.45 mm/minCrucible rotation rate: 5 RPMMagnetic field type: horizontalMagnetic field strength: 3000 GGas flow: no,Temperature field (K): top left: parallel cross section;bottom left: perpendicular cross section;

51、top right: top view,Case 1,2. Cz Si growth under a TMF (cont’d),FEMAGSoft © 2013,,Case 1,Crystal diameter: 400 mmPulling rate: 0.45 mm/minCrucible rotation rate: 5 RPMMagnetic field type: horizontalMagnetic fiel

52、d strength: 3000 GGas flow: no,Velocity field (m.s-1): top left: parallel cross section;bottom left: perpendicular cross section,2. Cz Si growth under a TMF (cont’d),FEMAGSoft © 2013,,,Case 1,Crystal diameter: 400

53、 mmPulling rate: 0.45 mm/minCrucible rotation rate: 5 RPMMagnetic field type: horizontalMagnetic field strength: 3000 GGas flow: no,Oxygen concentration (m-3): top left: parallel cross section in melt;bottom left:

54、perpendicular cross section in melt; top right: crystal,2. Cz Si growth under a TMF (cont’d),FEMAGSoft © 2013,,Case 2,Crystal diameter: 300 mmPulling rate: 0.45 mm/minCrucible rotation rate: 5 RPMMagnetic field t

55、ype: horizontalMagnetic field strength: 3000 GGas flow: no,,,Temperature field (K): top left: parallel cross section;bottom left: perpendicular cross section,,2. Cz Si growth under a TMF (cont’d),FEMAGSoft © 2013

56、,,Crystal diameter: 400 mmPulling rate: 0.45 mm/minCrucible rotation rate: 5 RPMMagnetic field type: horizontalMagnetic field strength: 3000 GGas flow: no,Temperature field (K): top left: parallel cross section;bot

57、tom left: perpendicular cross section,Case 1,2. Cz Si growth under a TMF (cont’d),FEMAGSoft © 2013,,Velocity field (m.s-1): top left: parallel cross section;bottom left: perpendicular cross section,Case 2,Crystal d

58、iameter: 300 mmPulling rate: 0.45 mm/minCrucible rotation rate: 5 RPMMagnetic field type: horizontalMagnetic field strength: 3000 GGas flow: no,2. Cz Si growth under a TMF (cont’d),FEMAGSoft © 2013,,Case 1,Crys

59、tal diameter: 400 mmPulling rate: 0.45 mm/minCrucible rotation rate: 5 RPMMagnetic field type: horizontalMagnetic field strength: 3000 GGas flow: no,Velocity field (m.s-1): top left: parallel cross section;bottom l

60、eft: perpendicular cross section,2. Cz Si growth under a TMF (cont’d),FEMAGSoft © 2013,,Oxygen concentration (m-3): top left: parallel cross section in melt;bottom left: perpendicular cross section in melt; top rig

61、ht: crystal,Case 2,,Crystal diameter: 300 mmPulling rate: 0.45 mm/minCrucible rotation rate: 5 RPMMagnetic field type: horizontalMagnetic field strength: 3000 GGas flow: no,2. Cz Si growth under a TMF (cont’d),FEMAG

62、Soft © 2013,,,Case 1,Crystal diameter: 400 mmPulling rate: 0.45 mm/minCrucible rotation rate: 5 RPMMagnetic field type: horizontalMagnetic field strength: 3000 GGas flow: no,Oxygen concentration (m-3): top left

63、: parallel cross section in melt;bottom left: perpendicular cross section in melt; top right: crystal,2. Cz Si growth under a TMF (cont’d),FEMAGSoft © 2013,,,Crystal diameter: 400 mmPulling rate: 0.45 mm/minCruci

64、ble rotation rate: 5 RPMMagnetic field type: horizontalMagnetic field strength: 5000 GGas flow: no,Temperature field (K): top left: parallel cross section;bottom left: perpendicular cross section; top right: top view

65、,Case 3,2. Cz Si growth under a TMF (cont’d),FEMAGSoft © 2013,,Crystal diameter: 400 mmPulling rate: 0.45 mm/minCrucible rotation rate: 5 RPMMagnetic field type: horizontalMagnetic field strength: 3000 GGas flo

66、w: no,Temperature field (K): top left: parallel cross section;bottom left: perpendicular cross section; top right: top view,Case 1,2. Cz Si growth under a TMF (cont’d),FEMAGSoft © 2013,,,Velocity field (m.s-1): top

67、 left: parallel cross section;bottom left: perpendicular cross section,Case 3,Crystal diameter: 400 mmPulling rate: 0.45 mm/minCrucible rotation rate: 5 RPMMagnetic field type: horizontalMagnetic field strength: 500

68、0 GGas flow: no,2. Cz Si growth under a TMF (cont’d),FEMAGSoft © 2013,,Case 1,Crystal diameter: 400 mmPulling rate: 0.45 mm/minCrucible rotation rate: 5 RPMMagnetic field type: horizontalMagnetic field strength

69、: 3000 GGas flow: no,Velocity field (m.s-1): top left: parallel cross section;bottom left: perpendicular cross section,2. Cz Si growth under a TMF (cont’d),FEMAGSoft © 2013,,,Oxygen concentration (m-3): top left: