How to speed up an ECLIPSE run: ECLIPSE Convergence Version 2 December 2009 1. Introduction .................................................................................................................................................................... 2 A case study.................................................................................................................................................................... 2 2. Background to timesteps and iterations .......................................................................................................................... 6 2.1 Reports of timesteps and iterations ........................................................................................................................... 7 3. Report Steps ................................................................................................................................................................... 8 4. Timesteps ....................................................................................................................................................................... 9 4.1 Maximum Timestep .................................................................................................................................................. 9 4.2 Maximum Timestep after a well change ................................................................................................................. 10 Summary: Using TUNING ........................................................................................................................................... 11 Rule 1 ....................................................................................................................................................................... 11 Rule 2 ....................................................................................................................................................................... 11 Rule 3 ....................................................................................................................................................................... 11 5. Non-Linear Convergence Criteria ................................................................................................................................ 12 5.1 A simple 1-variable example .................................................................................................................................. 13 5.2 Definition of convergence ...................................................................................................................................... 16 5.3 Definition of non-linear iterations .......................................................................................................................... 17 5.4 Two-variable example ............................................................................................................................................ 17 5.5 Definition of convergence ...................................................................................................................................... 21 5.6 Tracking the source of the problem ........................................................................................................................ 23 5.6.1 ECLIPSE 100 .................................................................................................................................................. 23 5.6.2 ECLIPSE 300 .................................................................................................................................................. 27 5.6.3 Discussion on Non-Linearity ........................................................................................................................... 30 5.6.4 Non-Linear Divergence ................................................................................................................................... 30 5.6.5 Summary: How to identify problem cells ........................................................................................................ 35 6. Improving the data ........................................................................................................................................................ 36 6.1 Messages,Comments, Warnings ............................................................................................................................. 36 6.2 Problems, Errors, Bugs and NaNs .......................................................................................................................... 37 6.3 What to look for in different Sections .................................................................................................................... 38 6.3.1 GRID data ....................................................................................................................................................... 38 6.3.2 SCAL data ....................................................................................................................................................... 40 6.3.3 PVT data.......................................................................................................................................................... 45 6.3.4 SOLUTION data ............................................................................................................................................. 48 6.3.5 SUMMARY section ........................................................................................................................................ 50 6.3.6 SCHEDULE data ............................................................................................................................................ 50 7. Summary: Convergence ................................................................................................................................................ 52 7.1 Convergence problem Checklist ............................................................................................................................. 53 1. Introduction This document describes how to improve the convergence and the speed of the ECLIPSE simulators by simple changes to the data file. The data chosen to model the reservoir may be such that the simulator can only solve the model by taking extremely short timesteps or excessive amounts of cpu time. Making small changes to the model or adding new options or new keywords can sometimes lead to dramatic improvements in the speed of the simulation without changing the results to engineering accuracy. We will show what changes can speed up the run in this way. If the data cannot be changed, convergence criteria can sometimes be chosen to improve simulator performance. You do not normally need to change convergence criteria, and we advise you not to do so unless necessary; internal parameters in ECLIPSE will usually give you a stable and robust solution for your model, in a reasonable cpu time, although of course a big complex model will take longer to run than a small simple model. In some cases, however, the results of a run can depend on the convergence criteria and how they are applied, and there is a trade-off between accuracy and speed. We will show how to detect these cases and how to choose suitable convergence criteria. Recognizing and correcting the cause of convergence problems is an important part of simulation. We will look at what you can do if a model fails to converge, or when it takes a long time to converge. We will show that: 1. Failure to converge is not (usually) due to a bug in ECLIPSE. It is usually related to the data. 2. Failure to converge Linear iterations is not always something to worry about. 3. Failure to converge a minimum timestep can affect the validity of the results. ECLIPSE can output a lot of messages, warnings, etc. You need to know what should you do about them, which ones you can safely ignore, and which ones need some action from you. You can also ask ECLIPSE to produce reports showing how both linear and non-linear iterations are proceeding and the methods by which time steps are selected. We will look at these reports. Most of the advice and suggestions apply equally to all the simulators. When there are differences in the detailed treatment between ECLIPSE Black Oil and ECLIPSE Compositional we will highlight these differences and explain how the data for each simulator can be tuned to improve performance. A case study We received an ECLIPSE model from a client. The model was “chopping” timesteps and the model was taking many hours to run, which was unexpected given the size of the model. We were asked if we could speed up the runs. We looked at the evidence. The main suspect was one of the production wells, so we looked at the data for that well in detail: -- WELL SPECIFICATION DATA -WELSPECS -- ---- WELL NAME GROUP LOCATION BHP PI NAME I J DEPTH DEFN PRODUCER G 10 10 8400 'OIL' / / -- COMPLETION SPECIFICATION DATA -COMPDAT --- WELL LOCATION OPEN/ SAT CONN -- NAME I J K1 K2 SHUT TAB FACT -PRODUCER 10 10 3 3 OPEN 0 -1 0.5 / PRODUCER 9 10 3 3 OPEN 0 -1 0.5 / PRODUCER 8 10 3 3 OPEN 0 -1 0.5 / PRODUCER 7 10 3 3 OPEN 0 -1 0.5 / PRODUCER 6 10 3 3 OPEN 0 -1 0.5 / PRODUCER 5 10 3 3 OPEN 0 -1 0.5 / PRODUCER 4 10 3 3 OPEN 0 -1 0.5 / PRODUCER 3 10 3 3 OPEN 0 -1 0.5 / PRODUCER 2 10 3 3 OPEN 0 -1 0.5 / / Both the well specification data and the completion data looked OK. The well production controls also looked OK. We looked at the “photographic” evidence. From the COMPDAT information we can see that the well is a horizontal well completed in layer 3. From the geometry information we knew that layer 3 was a dipping layer, so we expected to see a picture similar to What we saw was very different: Investigations showed that the run was a Restart run. The run had been part of in investigation as to whether the proposed well should be completed in layer 2 or layer 3. There had been some confusion over the names of the runs, and the base run had the well completed in layer 2 but the restart had the well completed in layer 3. ECLIPSE read the new COMPDAT data and concluded that it represented a workover, so the well now had completions in both layers. As there was poor communication between layer 2 and layer 3 and different pressures history in each, this created severe cross-flow. The cross-flow was the cause of the convergence issues and the poor performance. This short case study demonstrates that before we can solve a problem we need to know where the problem is. So what made us suspect this well? If a run is misbehaving, how do we find clues that will point us towards the source of the problem? More generally, how do we speed up an ECLIPSE run? What runs can be speeded-up? In fact all runs can be speeded up. There are three general solutions: 1. Hardware solutions - Get a faster machine. A new computer may be twice as fast as the equivalent computer from two years ago. - Make sure no other jobs are running on your computer. If you are running on a cluster, check that there is no contention for memory or other resources on the nodes on which your job is running. - Run in parallel. This will not speed up your run if there are convergence issues in serial mode. If your model is running smoothly but slowly on one processor, try running in parallel on two. If it speeds up, try running on four processors. 2. Reduce the size of the model - Can you reduce the number of grid blocks without affecting the quality of the solution? For instance can you replace a large number of waterfilled blocks with an analytic aquifer? - If you are running in compositional mode, can you model your fluid with fewer pseudo-components? 3. Identify data issues - Can you gain any speed by changing the time-stepping? - Can you identify the cause of any convergence problems? In the following chapter we will explain: • • • • How to get more information What the details of the output mean What to look for in the data How to fix some of the problems 2. Background to timesteps and iterations An ECLIPSE simulation A report step A time step A non-linear iteration is is is is made made made made up up up up of of of of one one one one or or or or more more more more report steps. time steps. non-linear iterations. linear iterations. When you set up the data model to run ECLIPSE, you are asked to specify the report steps that you want. You therefore have direct control of the number of report steps and the time gap between them. Reducing the number of report steps can sometimes reduce the cpu time. ECLIPSE has default values that control how many timesteps will be used to reach the next report that you have asked for. These default values will work well in most cases, but there are times when you may need to adjust some of the defaults so that the simulation takes fewer timesteps. In many (but not all) cases, fewer time steps will lead to less cpu time. In some difficult cases, reducing the maximum timestep may speed up the run. Different default values also control how many non-linear iterations will be used to solve each timestep. These values should normally be left unchanged. In a few cases, adjustments to the convergence criteria can improve the performance of the simulator. In most cases however the greatest improvements in performance are obtained by identifying the cause of the non-linear problem and changing the data model to reduce the non-linearity. The major part of this document will explain how to avoid problems of this type, and how to find and fix the problems if they do occur. By the time problems occur in the linear iterations, it is usually too late to fix them by adjusting the linear convergence control. Some controls can be changed in extreme circumstances, but the best advice is to avoid such problems by controlling the timestepping and the non-linear iterations. 2.1 Reports of timesteps and iterations The number of report steps, time steps and non-linear iterations can be found in both the PRT file and the LOG file - the shorter form of output which appears on the screen (for interactive runs) or is sent to a log file (for batch or background runs). On Linux (or Unix) systems, you can use the grep command to find all the necessary information; on PCs you may need to use your favourite editor and find the information one line at a time. Details are reported in a different way in ECLIPSE 100 and ECLIPSE 300. ECLIPSE 100 For a Linux system and an ECLIPSE 100 log file called BASE.LOG, the command grep TSTEP BASE.LOG >BASE.STEPS100 will create a file BASE.STEPS100 that will contain one line for each timestep. Each line will be of the form: STEP 15 TIME= 400.00 DAYS ( +30.0 DAYS REPT 3 ITS) (4-FEB-2003) "STEP 15" means this is the 15th timestep. "TIME= 400.00 DAYS" means there have been 400 simulated days since the beginning of the simulation. "+30.0 DAYS" shows that the latest timestep was of 30 days. "REPT" is a mnemonic explaining why 30 days were chosen. "REPT" means that a report step has been reached. "3 ITS" mean 3 non-linear iterations were needed to solve the 30 day timestep. "(4-FEB-2003)" is the current simulation date. For a number of timesteps, the BASE.STEPS100 file will look like this: STEP STEP STEP STEP STEP STEP STEP STEP STEP STEP STEP STEP 1 2 3 4 5 6 7 8 9 10 11 12 TIME= TIME= TIME= TIME= TIME= TIME= TIME= TIME= TIME= TIME= TIME= TIME= 2.40 HOURS ( +2.40 HOURS 9.60 HOURS ( +7.20 HOURS 1.30 DAYS (+21.60 HOURS 4.00 DAYS ( +2.7 DAYS 12.10 DAYS ( +8.1 DAYS 29.22 DAYS ( +17.1 DAYS 49.63 DAYS ( +20.4 DAYS 81.72 DAYS ( +32.1 DAYS 130.36 DAYS ( +48.6 DAYS 205.31 DAYS ( +75.0 DAYS 285.16 DAYS ( +79.8 DAYS 365.00 DAYS ( +79.8 DAYS INIT MAXF MAXF MAXF MAXF TRNC TRNC TRNC TRNC TRNC HALF REPT 4 2 2 3 3 4 4 4 5 3 3 3 ITS) ITS) ITS) ITS) ITS) ITS) ITS) ITS) ITS) ITS) ITS) ITS) (19-OCT-1982) (19-OCT-1982) (20-OCT-1982) (23-OCT-1982) (31-OCT-1982) (17-NOV-1982) (7-DEC-1982) (8-JAN-1983) (26-FEB-1983) (12-MAY-1983) (31-JLY-1983) (19-OCT-1983) ECLIPSE 300 For a Linux system and an ECLIPSE 300 log file called BASE.LOG, the command grep ";" BASE.LOG >BASE.STEPS300 will create a file BASE.STEPS300 that will contain one line for each timestep of the form: Rep ; 400.0 30.0 8.7838 .19498 1.4E05 32884. 1.2E06 4843.6 .00000 1.3E06 3 "Rep" is the mnemonic that shows that a report step has been reached, "30.0" shows that the latest timestep was of 30 days, the next 8 numbers show the GOR, water cut, oil/gas/water production rates, average field pressure, and gas and water injection rates, "3" at the end of the line means that 3 non-linear iterations were needed to solve the 30-day timestep. If the AIM option is used then the line will have and extra number at the end, Rep ; 400.0 30.0 8.7838 .19498 1.4E05 32884. 1.2E06 4843.6 .00000 1.3E06 3 2% 2% shows the percentage of the cells that was solved fully implicitly. For a number of timesteps, the BASE.STEPS300 file will look like this: Init MIF MIF MIF Rep Rep Rep Rep Rep Rep Rep Rep Rep Rep Rep Rep ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; 1.0 3.0 7.0 15.0 31.0 59.0 90.0 120.0 151.0 181.0 212.0 243.0 273.0 304.0 334.0 365.0 1.0 2.0 4.0 8.0 16.0 28.0 31.0 30.0 31.0 30.0 31.0 31.0 30.0 31.0 30.0 31.0 484.94 484.94 484.94 484.94 484.94 484.94 484.94 484.94 484.94 484.94 484.94 484.94 484.94 484.94 484.94 484.94 .01293 .01459 .01857 .02635 .04014 .05838 .07354 .08489 .09449 .10268 .11039 .11760 .12423 .13075 .13683 .14290 2.2E06 1.9E06 1.5E06 1.0E06 6.6E05 4.3E05 3.3E05 2.8E05 2.4E05 2.2E05 2.0E05 1.9E05 1.7E05 1.6E05 1.5E05 1.5E05 29221. 28463. 28320. 28039. 27512. 26824. 26277. 25872. 25539. 25271. 25037. 24834. 24658. 24494. 24350. 24213. 1.1E09 9.3E08 7.3E08 5.0E08 3.2E08 2.1E08 1.6E08 1.4E08 1.2E08 1.1E08 9.8E07 9.0E07 8.4E07 7.9E07 7.4E07 7.0E07 224.40 224.36 224.29 224.20 224.08 223.94 223.81 223.70 223.59 223.50 223.40 223.32 223.23 223.15 223.08 223.00 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 3 3 2 3 3 3 4 3 3 3 3 3 3 3 3 3 3% 3% 3% 3% 3% 3% 3% 3% 3% 3% 3% 3% 3% 3% 3% 3% Timestep Reason Some of the more common reasons for timestep selection are: Mnemonic E100 E300 Explanation INIT MAXF REPT HREP CHOP DIFF TRNC first timestep maximum increase factor report step half step to report timestep chopped follows CHOP TTE limit Solution Change Throughput Limit Init MIF Rep Hrep Redu TTE SCT TPT A more complete list is available in the ECLIPSE Technical Description manual, in the section on Convergence. 3. Report Steps The number of report steps and the time between report steps will depend on the type of model that you are simulating: For a prediction or forecasting run lasting for instance 30 years you may for instance ask for monthly reports for the first year, quarterly reports for the next 5 years, and yearly reports for the rest of the simulation. For history matching you may ask for weekly reports for the first year and for monthly reports for remainder of the history match, to test the validity of your model on a finer time scale. For slim-tube experiments, reporting intervals are likely to be minutes and hours Computing time may be reduced by changing the requested reports if the following three conditions apply: 1. You are asking for more reports than you actually need 2. Each report step is being reached in just a single timestep 3. Each time step is being solved in a small number of non-linear iterations. By small we ideally mean one, but savings could also be made if the timestep is taking 2 or 3 iterations. If you have created a file BASE.STEPS100 and find that it mainly contains lines of the form: STEP STEP STEP STEP 40 41 42 43 TIME= TIME= TIME= TIME= 400.00 410.00 420.00 430.00 DAYS DAYS DAYS DAYS ( ( ( ( +10.0 +10.0 +10.0 +10.0 DAYS DAYS DAYS DAYS REPT REPT REPT REPT 1 1 1 1 ITS) ITS) ITS) ITS) (4-FEB-2003) (14-FEB-2003) (24-FEB-2003) (6-MAR-2003) so that each timestep is a report step, and each timestep is solved in one iteration, then the run may go 2 or 3 times faster if you allow the simulator to produce report steps once a month instead of once every 10 days. The example above was for an ECLIPSE 100 run. 300 or 500 run with: The same applies for an ECLIPSE Rep Rep Rep Rep 1.2E06 1.2E06 1.2E06 1.2E06 ; ; ; ; 400.0 410.0 420.0 430.0 10.0 10.0 10.0 10.0 8.7838 8.7838 8.7838 8.7838 .19498 .19498 .19498 .19498 1.4E05 1.4E05 1.4E05 1.4E05 32884. 32884. 32884. 32884. 4843.6 4843.5 4843.4 4843.3 .00000 .00000 .00000 .00000 1.3E06 1.3E06 1.3E06 1.3E06 1 1 1 1 4. Timesteps Having reached one report step, ECLIPSE will decide what timestep to take next according to: how easy or difficult the previous timestep was how the convergence of the previous time step compared with the convergence targets how much simulation time is left until the next report step whether either you or the default limits forced it to take a particular time step You can set the timestep limits by using the TUNING keyword in any of the ECLIPSE simulators, or by using the TSCRIT keyword in ECLIPSE 300 or 500. 4.1 Maximum Timestep The first example of how a run can be speeded-up is similar to the report step example already given. If you have created a file BASE.STEPS100 and find that it mainly contains lines of the form: STEP STEP STEP STEP 40 41 42 43 TIME= TIME= TIME= TIME= 400.00 410.00 420.00 430.00 DAYS DAYS DAYS DAYS ( ( ( ( +10.0 +10.0 +10.0 +10.0 DAYS DAYS DAYS DAYS MAXS MAXS MAXS MAXS 1 1 1 1 ITS) ITS) ITS) ITS) (4-FEB-2003) (14-FEB-2003) (24-FEB-2003) (6-MAR-2003) The mnemonic MAXS means that this timestep is the maximum allowed according to the TUNING keyword. Again we see that each timestep is a report step, and each timestep is solved in one iteration, and the run may go 2 or 3 times faster if you allow the simulator to produce report steps once a month instead of once every 10 days. For ECLIPSE 300 the mnemonic is "Max", and you can set the maximum timestep either using the TUNING or the TSCRIT keyword. The syntax to set the maximum timestep to 30 days using the TUNING keyword is: TUNING 1* 30 / / / The TUNING keyword is available in both ECLIPSE 100 and in ECLIPSE 300. ECLIPSE 300 also has the TSCRIT and CVCRIT keywords; TSCRIT is timestep control and CVCRIT is convergence control. Approximately, TSCRIT + CVCRIT = TUNING To set the maximum timestep to 30 days using TSCRIT: TSCRIT 2* 30 / Note that the default maximum timestep in ECLIPSE 100 is 365 days and in ECLIPSE 300 is 50 days. The limit will only be 10 days if it has been set to 10 in a TUNING or TSCRIT keyword. The maximum timestep should be compatible with the report step interval. If you are asking for a report at the beginning of each month and have a timestep limit of 30 days, then those months with 31 days will need at least two timesteps. Increasing the timestep limit to 31 days would allow every month to be covered in one timestep if there were no other limitations. If you are asking for reports every 3 months, the maximum timestep should be at least 92 days. You may sometimes help the non-linear solver by setting a lower value for the maximum timestep size. If the simulator is frequently chopping timesteps and you cannot find the cause as described in section 5 below, then reducing the maximum timestep can speed up the run. For instance if the simulator converges all steps less than 20 days but chops steps over 20 days, then a 20 day limit will reduce computing time. You can always allow the timesteps to increase later in the simulation, after a difficult modelling problem has been overcome. You should however always first try to fix the cause of the problem before fixing the symptom. 4.2 Maximum Timestep after a well change A similar inefficiency can arise with the maximum step after a well change. There is a default limit of how big a timestep can be taken immediately after any well keyword has been used. The default value for ECLIPSE 300 is 20 days, and the default for ECLIPSE 100 is “unlimited”. If you are history matching and setting new well rates at the beginning of each month, then ECLIPSE 300 will be unable to take any timestep greater than 20 days. Rather than simulate January for example with a 20-day timestep followed by an 11-day timestep, ECLIPSE will even-out the timesteps by having two 15.5-day timesteps. The mnemonic for this choice of timestep is 'HALF' in ECLIPSE 100 and 'HRep' in ECLIPSE 300. The timestep summary will then be: STEP STEP STEP STEP STEP STEP STEP 1 2 3 4 5 6 7 TIME= TIME= TIME= TIME= TIME= TIME= TIME= 15.50 15.50 14.00 14.00 15.50 15.50 15.00 DAYS DAYS DAYS DAYS DAYS DAYS DAYS ( ( ( ( ( ( ( +15.5 +15.5 +14.0 +14.0 +15.5 +15.5 +15.0 DAYS DAYS DAYS DAYS DAYS DAYS DAYS HALF REPT HALF REPT HALF REPT HALF 1 1 1 1 1 1 1 ITS) ITS) ITS) ITS) ITS) ITS) ITS) (15-JAN-2002) ( 1-FEB-2002) (14-FEB-2002) ( 1-MAR-2002) (15-MAR-2002) ( 1-APR-2002) (15-APR-2002) Increasing the maximum step after a well change to be more than 31 days could improve the timestepping to be: STEP STEP STEP 1 2 3 TIME= TIME= TIME= 31.00 DAYS ( +31.0 DAYS REPT 28.00 DAYS ( +28.0 DAYS REPT 31.00 DAYS ( +31.0 DAYS REPT 1 ITS) 1 ITS) 1 ITS) ( 1-FEB-2002) ( 1-MAR-2002) ( 1-APR-2002) You may want to reduce the maximum timestep after a well change. If for instance you expect a major change in reservoir behaviour every time you change the well controls, then you can save cpu time by setting this value to for instance 1 day, using TUNING 9* 1 / / / or TSCRIT 10* 1 / Summary: Using TUNING Rule 1 Don‟t use TUNING (or TSCRIT or CVCRIT) Rule 1a: Remove these keywords if they are already there. Rule 2 OK, if you have to: You can change the initial, minimum, or maximum timestep If you change the minimum, remember to also change the minimum choppable step If you change the maximum, make it consistent with report steps (for instance have a maximum of 31 days not a maximum of 30 days) Rule 3 Convergence controls: remember rule 1 Look for most common reason for the selection of the timestep size Look for most common reason for a chopped timestep If you have poor convergence, it is sometimes better to tighten convergence controls than to weaken them. This is discussed below. 5. Non-Linear Convergence Criteria The equations that the simulators are trying to solve are non-linear. By nonlinear we mean that for instance doubling the tubing-head pressure of a water injector will not usually double the amount of water injected, and doubling the oil saturation in a grid block will not usually double the oil mobility in that grid block. The simulators use an iterative process based on Newton's method to solve these non-linear equations: 1. 2. 3. 4. We linearise the equations We solve the linear equations We check if this linear solution gives us a good enough non-linear solution. If it does then we move to the next timestep. If not we calculate the change needed to improve the solution, then go back to step 1. Advance Timestep Linearize the Equations Iterate to solve the linear equations Plug the linear solution into the non-linear equation “Non-linear iteration” No Is the solution good? Yes 5.1 A simple 1-variable example We can demonstrate the method for the simple case of a function of one variable only. Suppose F is a continuous monotonic function of x, and that F(x)=0 for some value of x. We want to find this value of x. We start by guessing a value for x. Let‟s call this first guess x0. We calculate F(x0) to see if F(x0)=0. In this case F(x0) is not zero so x0 is not the answer. Our next step is to “linearise” the problem by calculating the gradient of f(x) at x=x0. We cannot easily solve F(x)=0 so we solve an approximation to F(x), and the approximation here is a straight line. We solve to give us the value of x at which the straight line of the gradient cuts the x axis. This value is x1, and is an improvement on our initial guess of x0. We now calculate F(x1) to see if F(x1)=0. F(x1) is not zero, so x1 is not the answer we are looking for. Our next step is to “linearise” the problem again by calculating the gradient of f(x) at x=x1. We calculate the gradient of F at x1, solve the linear equation (the gradient equation) and obtain a new solution x2. We then check F(x2)to see if it is zero. As F(x2) is not zero, we continue our iterations and calculate the gradient of F at x2. This leads us to our third estimate of the solution, x3. Although F(x3) is not exactly zero, we decide that it is close enough for engineering accuracy, and accept x3 as our solution. 5.2 Definition of convergence The convergence criteria are different in ECLIPSE 100 and ECLIPSE 300. In ECLIPSE 100, the convergence criterion is as defined above. The absolute value of F has to be less than a defined limit. In ECLIPSE 300, the calculation of F is expensive as it involves flash calculations. For this reason, the criterion for convergence is that the change in x since the last iteration should be small. In the case of Dual Porosity ECLIPSE 100 runs, the TUNINGDP keyword will accept iterations as converged if either the residual or the solution change is small enough. The TUNINGDP keyword is sometimes helpful for high throughput cases. 5.3 Definition of non-linear iterations Each time the simulator goes through steps 1 to 3, it performs one non-linear iteration. The total number of times the simulator goes through these steps for each timestep is the number of non-linear iterations for that timestep. In the example given, there are 4 iterations, a first guess plus 3 improved values. There is a limit to the number of non-linear iterations that the simulators will try before giving up and trying with a smaller timestep. This limit depends on the simulator and on the solution method, and is set using either the TUNING or the CVCRIT keyword. You can check on how well the model is converging to a solution by looking at the number of non-linear iterations for each timestep, as described in section 2.1: - 1 non-linear per timestep means the step was very easy to converge - 2 to 3 non-linears per timestep means the step was easy to converge - 4 to 9 non-linears per timestep shows an increasingly difficult problem - > 10 non-linears per timestep can mean a problem with the model 5.4 Two-variable example We can generalise the one-variable example above to two variables, and demonstrate how the iteration procedure works in ECLIPSE for two variables. Consider a hill-climbing problem. Mathematically we have a function F of two variables x and y, and we want to find the value of x and the value of y for which F is a maximum. We can think of this as trying to climb to the top of a smooth rounded hill when it is so foggy that we can only see our feet. Suppose the hill is described by the contour lines below: We start at an initial position x0. x0 As with the 1-D case, we first check whether x0 is the solution. Mathematically we want to check if two orthogonal derivatives are zero, and if not then we want to calculate a planar approximation to the surface at x0. On our hill we do this by taking one small step in any direction, returning to our starting point, then taking another small step at right angles to the first. If the slopes in both directions are zero then we are at the top. In this case we are not at the top, but our two small steps will give us enough information to work out the steepest slope at x0. This gives us a direction for our first large step, which takes us to x1. x0 x1 At x1 we perform the same iterations. First we check if we are close enough to the solution. If we are not, we linearise the local surface by taking two small steps at right angles to each other, then we take a large step in the steepest vertical direction. We continue doing this until we are close enough to the solution. x0 x2 1 x3 x1 In a real case in ECLIPSE, we are dealing with 100,000s or millions of variables, the non-linear surface is much more complex, and the solution is more difficult to find. Solving the linear equations is also time-consuming when there are millions of variables and millions of simultaneous equations. These linear equations now require their own linear iterations, and this is where most of the computing time is likely to be spent in a reservoir simulation. A detailed description of the solution of the linear equations is beyond the scope of this manual. The interaction between the linear and non-linear equations is shown in the following picture: The horizontal axis represents time; the vertical axis represents the solution variables, which in the case of ECLIPSE 300 are the pressures and saturations in every grid blocks. We start with the known solution at time=n and we want to find the solution at time=n+1. Our first estimate of the solution at n+1 is that it is the same as the solution at n. As before, our first step is to check whether our current estimate of the solution is a good enough approximation to the solution. We find out it isn‟t. Our non-linear iteration counter l=1. We linearize the non-linear equations and create a linear equation of the form Ax=b, where A I known as the Jacobian. We solve the linear problem by iteration. Schematically we have shown this as a red line with 5 linear iterations. We now have a second estimate of the solution to the non-linear equations, indicated by the position l=2. We check to see if this is a solution, then linearise the equations again, solve the linear equations (indicated by 4 blue linear iterations) and reach our third estimate at l=3. We continue in this way until we converge the non-linear problem. The values at the end of the last non-linear iteration are the values we use to start the next timestep, starting at time=n+1. 5.5 Definition of convergence Formally, we need a test to decide when we can stop the iterations and carry on to the next timestep, i.e. to decide when the solution that we have is 'good enough'. A test of convergence can be based on either a small residual or a small change in solution (or both). The residual can be thought of as a measure of how close we are to solving the non-linear problem and is used in ECLIPSE 100 for the primary test. ECLIPSE 300 uses the change of solution as its basic convergence test. 5.5.1 Solution Change Check The variables used by ECLIPSE 300 are pressure and molar densities. ECLIPSE 300 calculates the maximum change over all cells for these variables and compares these maximum changes to two limits, one for pressure (SCONVP) and one for saturation (SCONVS). SCONVP SCONVS Set by 1st data item in CVCRIT Set by 7th data item in CVCRIT default 0.1 atm default 0.01 The pressure change is compared directly with SCONVP, while an approximation using field-wide properties converts molar density changes to effective saturation changes. These effective saturation changes are compared with SCONVS. The maximum solutions change over all cells and all solution variables is converted to an "actual-over-target" (aot) which is a ratio of the current maximum divided by the appropriate limit. When the value of aot falls below 1, then the timestep is converged. For example, suppose that the largest absolute pressure change in any grid block during the last non-linear iteration was 0.35 atm, with a limit on pressure change SCONVP=0.1 atm. The ratio of the two number is 3.5, meaning that the pressure change was 2.5 time bigger than we would accept. Suppose the limit on effective saturation change was 0.01 and the actual largest absolute value in any cell was 0.02, which means the actual change was twice what we would accept. The value aot is the maximum of 2.5 and 2, so aot=2.5. If the 8th data item in the DEBUG3 keyword is set greater than 0, then ECLIPSE 300 will output a line of debug that shows the convergence behaviour. The UNIX grep command can be used to produce summaries of the output. e.g. grep ";|aot" ROOT.DBG NLStep= NLStep= NLStep= NLStep= Rep ; NLStep= NLStep= NLStep= NLStep= NLStep= NLStep= MIF ; 0 lin= 1 lin= 2 lin= 3 lin= 8901.0 0 lin= 1 lin= 2 lin= 3 lin= 4 lin= 5 lin= 8903.0 will produce something like: 23 aot= 97.21 Rmax= .7162E-01 Rsum= .9919E-05 egain=-.1000E+01 19 aot= 17.55 Rmax= .1762E+00 Rsum= .1329E-06 egain=-.1000E+01 21 aot= 2.94 Rmax= .1421E-01 Rsum= .6285E-06 egain=-.1000E+01 12 aot= .77 Rmax= .7252E-02 Rsum= .7476E-08 egain=-.1000E+01 1.00 8.7838 .19498 1.4E05 32884. 1.2E06 4843.6 .00000 1.3E06 4 2% 27 aot= 137.50 Rmax= .7587E-01 Rsum= .2011E-04 egain= .1805E+00 26 aot= 79.23 Rmax= .7589E-01 Rsum= .1743E-04 egain= .1676E+00 26 aot= 76.18 Rmax= .7279E-01 Rsum= .1676E-04 egain= .2621E+00 24 aot= 9.30 Rmax= .9062E-01 Rsum= .8301E-06 egain=-.1000E+01 24 aot= 9.00 Rmax= .8764E-01 Rsum= .8028E-06 egain=-.1000E+01 13 aot= .08 Rmax= .9509E-01 Rsum= .9853E-09 egain=-.1000E+01 2.00 8.7860 .19501 1.4E05 32880. 1.2E06 4844.2 .00000 1.3E06 6 2% This shows two time steps: The first 4 lines above show iterations 0 to 3. Non-linear iterations start with iteration 0, which is a first guess at the new solution. The aot for iteration 0 is 97.21 so the non-linears are not converged. Iteration 3 has an aot of 0.77, which is less than 1. The largest absolute change is now less than the target, so the iterations have converged. The simulator provides a timestep report which starts “Rep” meaning a report has been reached, and a “4” near the end of that line means that 4 non-linears were needed to converge. The second timestep needs a total of 6 iterations. In this step we see that the aot at the last non-linear iteration drops from 9.00 to 0.08. This is a sign of “quadratic convergence” which will occur when the simulator is very close to the solution and reduces the aot by orders by magnitude at each non-linear iteration. Details If you wish to tighten (reduce) the convergence criteria, there are minimum values below which SCONVP and SCONVS cannot be set. The pressure minimum is 0.01 and the effective saturation minimum 0.005. The pressure minimum can be ignored by setting the data item (1st in CVCRIT) negative. This behaviour is not documented, and may be changed in future releases. The convergence tolerances are relaxed slightly at each non-linear iteration, to make convergence easier as the number of iterations increases. The hope is that this will allow a difficult timestep to be completed without significantly effecting the results. The factor used is: Factor = 1.0 + “iteration number”/”maximum iterations” An additional check is made that the sum of the residuals, which is a measure of total material balance error, is not excessive. This is a safety measure that can prevent convergence but rarely gets invoked in production cases. Even if the solution change is too great and aot>1, the maximum residual may be small enough for the timestep to be converged. The criterion used is the variable SNLRMX set using the 11th data item in CVCRIT. The residual considered is the residual calculated using the previous iteration and is not the new residual resulting from the solution change. Because of this the test tends not to be effective very often. 5.5.2 Gain Option in ECLIPSE 300 “Gain” is an option, not used by default, which can significantly improve performance. The idea is to speed the code up by not taking an extra non-linear iteration if that iteration is likely to generate a very small solution change. We predict the behaviour of the next iteration based on the history of previous timesteps. We calculate a gain factor from previous timesteps and rather than using the aot we use: EAOT=MIN(2.0*AOT*EGAIN,AOT) Where EAOT - Effective aot AOT - Is the aot as calculated above EGAIN - A measure of the expected improvement at the next step The gain value is printed out with the non-linear debug (see the example above). The gain option can be switched on using the 68th switch of the OPTIONS3 keyword. The gain option can be controlled using the 19th switch of the OPTIONS3 keyword. If set to 1, the gain will not be used if any cell has changed state during the iteration. In this case you might expect the previous history to be less valid. The option should be used with care. We have seen examples where this option can lead to dramatic changes in iteration sequence. 5.6 Tracking the source of the problem A single cell can cause non-convergence. As we increase the number of cells in a simulation, we increase the odds that a cell will cause non-convergence. If the are only one or two cells in the reservoir model that are causing problems, we can identify these cells and check if there is any engineering or data reason which could explain why they are causing problems. For example the cells may be at or near well completions, in which case well control could be modified, or could be cells with very small pore volumes in which case the MINPV keyword could be used. 5.6.1 ECLIPSE 100 In ECLIPSE 100, you will get non-linear debug if you set “NEWTON=2” in RPTSCHED. This will produce output of the form: Looking in detail at an example of the residuals, IT= 0 CNV CELL MAT BAL OIL 1.00424( 28, 45, 3) 5.3D-03 WAT-0.00288( 9, 3, 3) -1.3D-07 GAS********( 5, 45, 1) -1.3D-02 DPRESS 0.00 0.00 0.00 DSWAT 0.00000 0.00000 0.00000 DSGAS 0.00000 0.00000 0.00000 LINIT= 5 NSCHP= 6 NCHOP= IT= 1 CNV CELL MAT BAL OIL-1.99144( 5, 45, 1) 4.3D-02 WAT-0.16316( 2, 4, 4) -3.7D-06 GAS********( 5, 45, 1) -3.1D-02 0 NSTAT1,2,3= 50 5400 DPRESS DSWAT DSGAS -24.89 0.00026 -0.20000 -14.02 0.00490 0.00000 -24.89 0.00026 -0.20000 0 NTRAN= 321 LINIT= 3 NSCHP= 195 NCHOP= IT= 2 CNV CELL MAT BAL OIL-0.62319( 5, 45, 1) 1.7D-02 WAT-0.04162( 28, 5, 3) -1.3D-04 GAS********( 5, 45, 1) -2.2D-02 0 NSTAT1,2,3= 50 5370 DPRESS DSWAT DSGAS -21.15 0.00081 -0.01843 -30.47 0.00139 0.04000 -21.15 0.00081 -0.01843 30 NTRAN= 30 LINIT= 3 NSCHP= 44 NCHOP= IT= 3 CNV CELL MAT BAL OIL-0.30993( 5, 45, 1) -8.9D-05 WAT-0.04591( 28, 4, 3) -4.0D-05 GAS********( 5, 45, 1) -1.4D-02 3 NSTAT1,2,3= 50 5367 DPRESS DSWAT DSGAS -26.44 0.00088 -0.21134 -19.80 0.00666 0.11687 -26.44 0.00088 -0.21134 33 NTRAN= 25 This shows the first 4 non-linear iterations (IT=0, IT=1, IT=2, IT=3) in a case that has convergence problems for the gas phase. The first line shows IT=0, the first iteration, and CNV etc are column headers for the next 3 lines. The columns are: CNV CELL MAT BAL DPRESS DSWAT DSGAS The The The The The The worst residual for the OIL, WATer and GAS phases cell that has the worst residual material balance for that cell, a measure of mass accuracy change in pressure in that cell since the last iteration change in water saturation since the last iteration change in gas saturation since the last iteration The residual for gas in all 4 iterations is shown as ******* which means that it is greater than the maximum printable value. It has a very high residual at each iteration for cell (5,45,1), so that is the cell that is causing problems. After each iteration report above, the line starting LINIT= provides more information on what is happening within the model. LINIT number of iterations required to solve the linearized equations. NSCHP number of saturation changes that were altered to suppress possible oscillations. NCHOP number of times the changes in P, Rs, or Rv were reduced to increase stability. NSTAT1,2,3 is the number of cells in solution state 1,2,3 Solution state 1 means no oil is present in the cell Solution state 2 means both oil and gas are present in the cell Solution state 3 means no gas is present in the cell NTRAN is the number of state transitions since the last non-linear iteration. Any non-zero value of NSCHP or NCHOP increases material balance errors for the subsequent non-linear iteration and therefore reduces the chances of convergence. Some saturation chops can be avoided by adjusting relative permeability curves in such a way that the critical saturation is not the same as the lowest saturation value in the table. For instance, instead of SWFN 0.2 0.3 0.4 0.5 0.6 0.8 0.9 1 0 0.07 0.15 0.24 0.33 0.65 0.83 1 7 4 3 2.5 2 1 0.5 0 / 0 0 0.07 0.15 0.24 0.33 0.65 0.83 1 7 1* 4 3 2.5 2 1 0.5 0 / try using SWFN 0.2 0.21 0.3 0.4 0.5 0.6 0.8 0.9 1 The new saturation value at 0.21 may help convergence. It will not affect the initial fluids-in-place but will unfortunately slightly reduce the water mobility for water saturations between 0.2 and 0.3. This may not be important to engineering accuracy. Look for oscillations in the CNV for a phase. If one iteration has a positive value, the next iteration has a negative value, then the next is positive, then negative, etc. then there is perhaps a non-linearity in the system. These are sometimes associated with sudden changes in the slope of the relative permeability curves. If you have access to the SCAL program, you can plot these slopes and look for discontinuities. If you have access to a spreadsheet program then you can numerically calculate and plot the slopes. Remember that ECLIPSE will use all the values of saturation and relative permeability that you give in the table without any smoothing. You should therefore try to avoid tables such as SWFN 0.2 0.21 0.3 0.301 0.398 0.4 0.401 0.402 0.5 0.6 0.8 0.9 1 0 0 0.07 0.07 0.14 0.15 0.17 0.19 0.24 0.33 0.65 0.83 1 7 1* 4 4 3 3 3 3 2.5 2 1 0.5 0 / The table above has saturation values that are too close to each other and the slopes of the relative permeabilities shows severe changes. You should also try to avoid tables such as SWFN 0.2 0.3 0.5 0.51 0.8 0.9 1 0 0. 0.01 0.60 0.68 0.83 1 7 4 2.5 2 1 0.5 0 / The table above has a very sudden krw change from krw=0.01 at Sw=0.5 to krw=0.60 at Sw=0.51 and will certainly cause convergence problems. TUNINGDP output The NEWTON switch in RPTSCHED will produce extra information in the case of TUNINGDP: PTRG is the target pressure change; the default is 1 psi in Field units. STRG is the target saturation change; the default is 0 if TUNINGDP is not used and 0.01 if it is used. MDDP is the maximum Pressure change for convergence. MDDS is the maximum Saturation change for convergence. If you use TUNINGDP, (i) ECLIPSE will solve the linear equations to a tighter tolerance (ii) the convergence is reached if either the residual or the solution change criteria is small enough. If you don‟t use TUNINGDP, only the residual is used to test for convergence. Output of reason for non-linear failure If you set DEBUG item 1 in ECLIPSE 100 to be > 1, ECLIPSE will output a number that shows the reason that a non-linear iteration has failed to converge. The 5 possible values are: 1 2 3 4 5 – Exceeds maximum tolerable pressure change DDPLIM Exceeds maximum tolerable saturation change DDSLIM Exceeds maximum tolerable material balance RSUM Exceeds tolerable for maximum norm RNMAX Uses but fails the TRGDDP/TRGDDS test DDPLIM and DDSLIM are described in record 3 of the TUNING keyword. A value of 3 means that the normalised residual is greater than the maximum allowed. That allowed maximum is a linear combination of items 3 and 7 of record 2 of the tuning parameter. The maximum residual also depends on which newton iteration we are on. A value of 4 means that the normalised residual is greater than the maximum allowed. That allowed maximum is a linear combination of items 2 and 6 of record 2 of the tuning parameter. The maximum residual also depends on which newton iteration we are on. A value of 5 means that the change in the solution is less than TRGDDP and TRGDDS. These are the minimum values of DDPLIM and DDSLIM until that point in the simulation. This is used in cases where there is high throughput. 5.6.2 ECLIPSE 300 In ECLIPSE 300, you can have either a visual display of the grid blocks causing convergence problems or you can look at numerical output in the same way as for ECLIPSE 100. Starting with the 2009.1 version, ECLIPSE 100 may also have the visual display option. Adding CONV to the RPTRST keyword will send two new outputs to the restart files. Each cell will have two new variables that will be used to count the number of times that cell has been one of the most difficult cells to converge. At the beginning of the simulation each cell will have its counter set to zero. At every timestep, the 10 most difficult cells will have their counter increased by 1. At the end of the run you can display the cells with the most problems. Visual option in RPTRST – ask for CONV Output of cells which are causing convergence problems. By default, CONV=10 is set so that the worst 10 cells will be output Output CONV_VBR: Worst cells based on volume balanced residual CONV_PRU: Worst cells based on pressure updates Alternatively, you will get non-linear debug if the 8th data item in the DEBUG3 keyword is set greater than 0. Some of the debug information has already been described in section 5.2 above. Typical output is of the form: Iteration 0 linears req 7 DX Pressure 0 -40.075969 25 32 4 F 1.469590 DX Comp 1 -0.000375 25 32 1 T 0.010000 DX Comp 2 -0.000980 25 32 1 T 0.010000 DX Comp 3 -0.025419 25 32 1 F 0.010000 DX Comp 4 -0.002318 25 32 1 T 0.010000 DX Comp 5 -0.000318 25 32 3 T 0.010000 DX Comp 6 0.009711 25 32 1 T 0.010000 DX Comp 7 0.008562 25 32 1 T 0.010000 DX Comp 8 0.004396 25 32 1 T 0.010000 DX Comp 9 0.001465 25 32 1 T 0.010000 DX Comp 10 0.000907 25 32 4 T 0.010000 NLStep= 0 lin= 7 aot= 27.27 Rmax=0.8514E+00 Rsum=0.2739E-03 egain=0.2624E-01 Iteration 1 linears req 5 DX Pressure 0 -0.563835 25 32 1 T 1.592056 DX Comp 1 0.000035 25 32 1 T 0.010833 DX Comp 2 0.000063 25 32 1 T 0.010833 DX Comp 3 0.001887 25 32 1 T 0.010833 DX Comp 4 0.000242 25 32 1 T 0.010833 DX Comp 5 0.000111 25 32 1 T 0.010833 DX Comp 6 0.000270 26 32 1 T 0.010833 DX Comp 7 0.000238 26 32 1 T 0.010833 DX Comp 8 0.000127 26 32 1 T 0.010833 DX Comp 9 0.000044 26 32 1 T 0.010833 DX Comp 10 -0.000317 25 32 2 T 0.010833 NLStep= 1 lin= 5 aot= 0.38 Rmax=0.1921E-01 Rsum=0.1448E-04 egain=0.3165E-01 Max changes:pres 40.6 25 32 4 temp 0.00 0 0 0 oil satn 0.516E-01 17 7 1 gas satn -0.178E-01 17 8 1 wat satn -0.939E-03 25 32 4 eng dens 0.00 0 0 0 Throughput ratio:avrg 0.404E-01 max 0.192 26 32 2 MIF ; 103.0 9.00 6.9094 .01725 8973.3 157.46 62000. 3531.4 0.0 60500. 2 2% This output shows two non-linear iterations leading to a timestep report. The first line: Iteration 0 linears req 7 tells us that the iteration 0 (the initial estimate) needed 7 linear iterations to solve the linear problem. The next 11 lines consist of one line for each of the solution variables showing the largest change in that variable in any cell during this iteration. The solution variables for each grid block are the pressure in the grid block and the molar density of each hydrocarbon component, and a water term. This model has 9 hydrocarbon components. Water is written as component 10. DX Pressure 0 DX Comp 1 DX Comp 2 DX Comp 3 DX Comp 4 DX Comp 5 DX Comp 6 DX Comp 7 DX Comp 8 DX Comp 9 DX Comp 10 -40.075969 -0.000375 -0.000980 -0.025419 -0.002318 -0.000318 0.009711 0.008562 0.004396 0.001465 0.000907 25 25 25 25 25 25 25 25 25 25 25 32 32 32 32 32 32 32 32 32 32 32 4 1 1 1 1 3 1 1 1 1 4 F T T F T T T T T T T 1.469590 0.010000 0.010000 0.010000 0.010000 0.010000 0.010000 0.010000 0.010000 0.010000 0.010000 In each of these lines, DX means the solution change. The first line is the Pressure change. The largest pressure change was an increase of 40.075969 psi in cell (25,32,4), which happen to contain an injecting completion. The “F” on that line means “False”, in that the pressure variable has not converged, as the pressure change is greater than 1.469590, which is the maximum pressure change, allowed for convergence for this iteration. The second DX line shows the largest change in the molar density, expressed as a saturation equivalent, for component 1. This increase of 0.000375 was in cell (25,32,1) and is less than the convergence maximum of 0.01, so that the component 1 variable is considered to be converged. In fact all the components have converged except for component 3. The non-linear iteration however has not converged since two of the variables (pressure and component 3) are not yet converged. The next line is a summary of the first iteration (iteration 0) and has already been discussed in section 5.5 above. NLStep= 0 lin= 7 aot= 27.27 Rmax=0.8514E+00 Rsum=0.2739E-03 egain=0.2624E-01 Nltep=0 states that this in non-linear step 0 Lin=7 says that 7 linear iterations were needed to solve it Rmax= 0.8514E+00 is the worst (maximum) residual at the beginning of this iteration Rsum= 0.2739E-03 is the sum of all the residuals at the beginning of this iteration egain=0.2624E-01 is the gain discussed in section 5.5.2 The next line Iteration 1 linears req 5 is the start of the report on the second non-linear iteration, Iteration 1, which needed 5 linear iterations. The next 11 lines are the solution changes. They are similar to the reported changes for the first non-linear iteration except that now both the pressure and component 3 have changed by less than the new convergence criteria. The nonlinear iterations have now converged. NLStep= 1 lin= 5 aot= 0.38 Rmax=0.1921E-01 Rsum=0.1448E-04 egain=0.3165E-01 is the report of the converged iteration showing that the maximum residual at the beginning of that iteration was down to 0.1921E-01 and the sum of residuals was down to Rsum=0.1448E-04. There then follows a report on the changes during that timestep. Max changes:pres 40.6 oil satn 0.516E-01 wat satn -0.939E-03 25 32 17 7 25 32 4 temp 0.00 1 gas satn -0.178E-01 4 eng dens 0.00 0 17 0 0 8 0 0 1 0 The maximum pressure and saturation changes are reported, as well as the cells in which this change occurred. In the case of thermal runs, the maximum temperature and energy density changes are also reported. The maximum throughput flowing through a cell throughput is too high problems, and the pore examined. Throughput ratio:avrg is reported next. Throughput is defined as the volume divided by the pore volume of the cell. If the (say higher than 0.5) it could cause convergence volume of the cell that has the high throughput should be 0.404E-01 max 0.192 26 32 2 The last of these lines is the report of production, etc. for the timestep. These report steps have been described in section 3. MIF ; 103.0 9.00 6.9094 .01725 8973.3 157.46 62000. 3531.4 0.0 60500. 2 2% Moving goalposts ECLIPSE 300 non-linear convergence criteria have 'moving goalposts'. The convergence tolerances are relaxed slightly at each non-linear iteration. Effectively, convergence becomes easier as the number of iterations increases. The idea is that if you have reached the maximum number of iterations (call that Nmax) and you are close (within a factor of two) to the convergence criteria, then you don't want to chop the timestep and waste all the work you have done so far. So the criteria relaxes by 1/Nmax every Newton iteration. An example is shown below. The maximum number of non-linear iterations is 12, the units are metric, and the first Newton iteration uses the default convergence criteria (0.1 atm pressure and 0.01 for component specific volume). For later Newton iterations, these criteria were relaxed by 8.33% (=1/12) with each Newton iteration. After 12 Newton iterations, the criteria are doubled to 0.2 atm pressure and 0.02 for component specific volume. Iteration 0 linears req 4 NTOTNL 4076 DX Pressure 0 2.803147 38 16 10 Ind:GLOB F 1.469590 DX Comp 1 0.000262 34 31 6 Ind:GLOB T 0.010000 DX Comp 2 0.000042 34 31 6 Ind:GLOB T 0.010000 .. NLStep= 0 lin= 4 aot= 1.91 Rmax=0.1149E-02 Rsum=0.8487E-05 egain=0.9926E+00 DCHOP2: 1 cells chopped, Try= 1 Iteration 1 linears req 7 NTOTNL 4077 DX Pressure 0 -2.799664 38 16 10 Ind:GLOB F 1.592056 DX Comp 1 -0.003126 38 16 10 Ind:GLOB T 0.010833 DX Comp 2 -0.000239 38 16 10 Ind:GLOB T 0.010833 .. NLStep= 1 lin= 7 aot= 1.91 Rmax=0.3573E-01 Rsum=0.1269E-05 egain=0.3326E+01 DCHOP2: 1 cells chopped, Try= 1 Iteration 2 linears req 4 NTOTNL 4078 DX Pressure 0 2.694699 38 16 10 Ind:GLOB F 1.714522 DX Comp 1 -0.048358 38 16 10 Ind:GLOB F 0.011667 DX Comp 2 -0.002969 38 16 10 Ind:GLOB T 0.011667 5.6.3 Discussion on Non-Linearity A few general problems are: 1. In explicit cells, the algorithm can try to extract more fluid from block than is present in the grid block. This happens because flows calculated from the solution of the previous timestep, and not from current timestep. The only thing that can be done is to reduce the a grid are the timestep. This will only be an issue if you using the IMPES or AIM solution methods. The default in ECLIPSE 100 is FULLIMP, but the default in ECLIPSE 300 is AIM unless you are using some options such as Radial, Dual Porosity or Thermal. 2. Flow reversals are a major non-linearity. 3. The group control algorithm sometimes wants to change well rates at every non-linear iteration to reflect the latest calculated conditions in the reservoir. Changing well rates every non-linear iteration can however lead to poor convergence. The simulators therefore by default only recalculate group control parameters for the first 4 non-linear iterations, then keep the group controls unchanged for the remaining non-linears. This value of 4 can be changed using NUPCOL, but we recommend leaving the value unchanged. 4. Non-monotonic VFP tables can cause convergence problems. VFP tables are always checked for monotonicity in ECLIPSE 300, and the check can be switched on in ECLIPSE 100 by using the EXTRAPMS keyword. 5.6.4 Non-Linear Divergence Earlier we showed a simple example of how the non-linear iteration converged for a simple function F(x). That function F(x) is “smooth” in the sense that there are no sudden changes in the values of the derivatives. What if the curve is not smooth? Consider the following function: We want to find the solution (marked by the word “Root”) but the function Y has a large slope change on either side of the solution. As usual, we start with a guess x=X0 and calculate the value of Y at x=X0. We find that X0 is not the solution, so we linearise Y by taking the slope and solving the linear equation of the slope to calculate a new value X1. We check X1 and find that X1 is not a solution either, so we have to continue the non-linear iterations. Note that X1 is further away from the correct solution at “Root” than X0, so that the non-linear iteration is moving us away from the solution. This next iteration takes us to a value X2, which is even further away from the correct solution. X3 is even further away, and the divergence continues ... If the data causes curves of this shape, then after a lot of computing time ECLIPSE will reach the limit on the number of non-linear iterations, chop the time step, and start iterating again. Unfortunately a smaller timestep may not solve the problem: a smaller timestep may mean that the initial guess is nearer the solution, but the slope change is still there and can still cause divergence. In some cases, the non-linear iterations can oscillate between two solutions: This oscillation or “flip-flop” can be seen in the debug output: LINIT=10 NSCHP= 6 NCHOP= IT= 6 CNV CELL MAT BAL OIL -0.03205( 10, 11, 3) 1.1D-06 WAT-0.00004( 10, 11, 3) -2.6D-08 GAS-0.07982( 13, 10, 3) -3.5D-06 LINIT= 6 NSCHP= 0 NCHOP= IT= 7 CNV CELL MAT BAL OIL -0.00004( 10, 18, 3) 5.0D-09 WAT 0.00000( 8, 17, 3) -1.0D-13 GAS 0.00354( 13, 10, 3) 2.0D-09 0 NSTAT1,2,3= 0 1793 1213 NTRAN= DPRESS DSWAT DSGAS -84.80 0.00002 -0.01805 -84.80 0.00002 -0.01805 -82.33 0.00006 -0.12626 0 NSTAT1,2,3= 0 1793 1213 NTRAN= DPRESS DSWAT DSGAS -84.65 0.00006 -0.00989 -86.77 0.00006 -0.04626 -81.34 0.00006 -0.11912 LINIT=10 NSCHP= 6 NCHOP= IT= 8 CNV CELL MAT BAL OIL -0.03205( 10, 11, 3) 1.1D-06 WAT-0.00004( 10, 11, 3) -2.6D-08 GAS-0.07982( 13, 10, 3) -3.5D-06 LINIT= 6 NSCHP= 0 NCHOP= IT= 9 CNV CELL MAT BAL OIL -0.00004( 10, 18, 3) 5.0D-09 WAT 0.00000( 8, 17, 3) -1.0D-13 GAS 0.00354( 13, 10, 3) 2.0D-09 0 NSTAT1,2,3= 0 1793 1213 NTRAN= DPRESS DSWAT DSGAS -84.80 0.00002 -0.01805 -84.80 0.00002 -0.01805 -82.33 0.00006 -0.12626 0 NSTAT1,2,3= 0 1793 1213 NTRAN= DPRESS DSWAT DSGAS -84.65 0.00006 -0.00989 -86.77 0.00006 -0.04626 -81.34 0.00006 -0.11912 3 0 3 0 The solution for non-linear iteration 6 is the same as that for non-linear iteration 8, and the solution for non-linear iteration 7 is the same as that for non-linear iteration 9. LINIT=10 NSCHP= 6 NCHOP= IT= 6 CNV CELL MAT BAL OIL -0.03205( 10, 11, 3) 1.1D-06 WAT-0.00004( 10, 11, 3) -2.6D-08 GAS-0.07982( 13, 10, 3) -3.5D-06 LINIT= 6 NSCHP= 0 NCHOP= IT= 7 CNV CELL MAT BAL OIL -0.00004( 10, 18, 3) 5.0D-09 WAT 0.00000( 8, 17, 3) -1.0D-13 GAS 0.00354( 13, 10, 3) 2.0D-09 0 NSTAT1,2,3= 0 1793 1213 NTRAN= DPRESS DSWAT DSGAS -84.80 0.00002 -0.01805 -84.80 0.00002 -0.01805 -82.33 0.00006 -0.12626 0 NSTAT1,2,3= 0 1793 1213 NTRAN= DPRESS DSWAT DSGAS -84.65 0.00006 -0.00989 -86.77 0.00006 -0.04626 -81.34 0.00006 -0.11912 LINIT=10 NSCHP= 6 NCHOP= IT= 8 CNV CELL MAT BAL OIL -0.03205( 10, 11, 3) 1.1D-06 WAT-0.00004( 10, 11, 3) -2.6D-08 GAS-0.07982( 13, 10, 3) -3.5D-06 LINIT= 6 NSCHP= 0 NCHOP= IT= 9 CNV CELL MAT BAL OIL -0.00004( 10, 18, 3) 5.0D-09 WAT 0.00000( 8, 17, 3) -1.0D-13 GAS 0.00354( 13, 10, 3) 2.0D-09 0 NSTAT1,2,3= 0 1793 1213 NTRAN= DPRESS DSWAT DSGAS -84.80 0.00002 -0.01805 -84.80 0.00002 -0.01805 -82.33 0.00006 -0.12626 0 NSTAT1,2,3= 0 1793 1213 NTRAN= DPRESS DSWAT DSGAS -84.65 0.00006 -0.00989 -86.77 0.00006 -0.04626 -81.34 0.00006 -0.11912 3 0 3 0 5.6.5 Summary: How to identify problem cells 1. Add „CONV‟ to RPTRST 2. To get more information for ECLIPSE 100, add „NEWTON=nn‟ to RPTSCHED 3. To get more information for ECLIPSE 300, DEBUG3 7* 1 / 6. Improving the data Make sure that the MESSAGES keyword has not been used to limit the number of messages (and comments, warnings, etc.) output by the program 6.1 Messages, Comments, Warnings First check all messages, comments and warnings. All of these start with an “@”, so you can find all these messages by looking for “@”, either in the PRT file or in the LOG file if you are running in batch mode. Alternatively ECLIPSE Office or PETREL RE will list all of these for you. Messages give information. You do not have to take any action about messages, but have a look at each message to make sure that the information is not unexpected. Two example are: @--MESSAGE AT TIME 0.0 DAYS (19-OCT-1982): @ CHECKING FOR LICENSES @--MESSAGE AT TIME 0.0 DAYS (19-OCT-1982): @ THE MEMORY REQUIRED TO PROCESS THE GRID SECTION IS @ 1526 KBYTES Comments also give information. You should read all comments to confirm that the model is what you want it to be. An example is: @--COMMENT @ AT TIME 0.0 DAYS (19-OCT-1982): NO NON-NEIGHBOUR CONNECTIONS FOUND If you did not expect any non-neighbor connections, then you do not have to do anything. If however you know that your reservoir has a lot of fault throws, then you may need to check your model. Warnings are things you should check. There may be data issues that need correcting. Common Warning (1) @--WARNING @ @ @ @ AT TIME 0.0 DAYS (19-OCT-1982): INCONSISTENT END POINTS IN SATURATION TABLE THE MAXIMUM GAS SATURATION (1.0000) PLUS THE CONNATE WATER SATURATION (0.1200) MUST NOT EXCEED 1.0 1 Probably safe, but your gas and water SCAL tables are inconsistent. Did you know that with your data the gas saturation will never be > 0.88? Is this what you wanted? Common Warning (2) @--WARNING @ @ AT TIME 0.0 DAYS (19-OCT-1982): THE BOTTOM HOLE PRESSURE LIMIT FOR WELL INJECTR1 HAS BEEN DEFAULTED. THE DEFAULT VALUE IS 100000. PSIA Probably safe as long as the injector never goes to BHP control, but setting a more reasonable limit would be safer. 6.2 Problems, Errors, Bugs and NaNs All of these also start with an “@”. Check and correct any problems and errors. If after you have fixed all problems and errors you still have either a bug or a NaN then contact your local support. Problems mean there is definitely an issue. Check and fix before continuing. @--PROBLEM @ @ @ @ @ AT TIME 0.0 DAYS (19-OCT-1982): OIL PRESSURE IN EQUILIBRATION CALCULATION BETWEEN DEPTHS 8400.0 FEET AND 8397.7 FEET HAS NOT CONVERGED. CONVERGENCE ERROR = 0.86508E+33 PSI . THE CURRENT NUMBER ( 100 ) OF DEPTH NODES IN THE EQUILIBRATION CALCULATION IS TOO SMALL. In this case, pointing to oil PVT table. Errors mean you cannot continue without fixing your data. @-@ @ @ ERROR AT TIME 0.0 DAYS (19-OCT-1982): ERROR IN PVTO TABLE NUMBER 1 NOT ENOUGH PRESSURE VALUES (ONLY 1) SPECIFIED FOR HIGHEST RS (= 1.270 MSCF/STB) This was the reason for the previous message, so always fix the first error message first. NaN means “Not a Number”, and is usually caused by dividing by zero or something similar. NaNs sometimes happens after an Error, so if you see “NaN” anywhere in any output, then: 1. Fix any Error or Problem first. 2. If the NaN is still there, contact your local support. @--WARNING @ @ @ @ @ @ @ AT TIME 0.0 DAYS (19-OCT-1982): GAS IS DENSER THAN OIL WITHIN THE RESERVOIR IN EQUILIBRATION REGION 1. CHECK SURFACE DENSITIES, FORMATION VOLUME FACTORS, RS, RV. THE PROBLEM OCCURS AT DEPTH 8725.0000 GAS PRESSURE = 0.4201E-01 OIL PRESSURE = NaN (RS MAY BE SPECIFIED IN THE WRONG UNITS) In this case, the NaN follows an error message. 6.3 What to look for in different Sections 6.3.1 GRID data MINPV Use MINPV to remove cells with small pore volumes. If a full-field model has convergence problems and does NOT have a minimum pore volume, adding MINPV is the most likely change that will improve performance. Take care to choose a suitable value of MINPV that does not significantly change the total pore volumes in each region, and that does not remove highpermeability streaks or thin shale layers within the reservoirs. If MINPV is being used to remove pinch-outs, also use PINCH to connect across those thin cells that represent pinch-outs. Within these limitations, using MINPV should both improve performance and give unchanged production results. LIMITS Look for the table headed “LIMITS” in the ECLIPSE 300 PRT file. This is a table of minimum and maximum values for each GRID array. Check the ratio of Maximum/Minimum values at the left of output, and investigate if any of these ratios is very large. Grid Array Maximum value Minimum value Non-zero minimum Ratio max/min ---------- ---------------------------- ---------------------------- ---------------------------- -------------TOPS 7571.0 ( 1 5 9) 7021.0 ( 1 1 1) 7021.0 ( 1 1 1) 1.078336 MIDS 7926.0 ( 1 5 9) 7146.0 ( 1 1 1) 7146.0 ( 1 1 1) 1.109152 DEPTH 7926.0 ( 1 5 9) 7146.0 ( 1 1 1) 7146.0 ( 1 1 1) 1.109152 DX 220.0 ( 1 1 1) 55.0 ( 1 6 2) 55.0 ( 1 6 2) 4.0 DY 220.0 ( 1 3 1) 55.0 ( 1 1 1) 55.0 ( 1 1 1) 4.0 DZ 710.0 ( 1 1 3) 250.0 ( 1 1 1) 250.0 ( 1 1 1) 2.84 DZNET 710.0 ( 1 1 3) 250.0 ( 1 1 1) 250.0 ( 1 1 1) 2.84 NTG 1.0 ( 1 1 1) 1.0 ( 1 1 1) 1.0 ( 1 1 1) 1.0 PORO 0.10 ( 1 1 1) 0.15E-01 ( 1 3 4) 0.15E-01( 1 3 4) 6.666667 PORV 459036.7 ( 1 1 3) 323.26 ( 1 1 4) 323.26 ( 1 1 4) 1420.0 PERMX 300.0 ( 1 1 4) 0.0 ( 1 1 1) 1.5 ( 1 1 6) 200.0 PERMY 300.0 ( 1 1 4) 0.0 ( 1 1 1) 1.5 ( 1 1 6) 200.0 PERMZ 7.5 ( 1 1 4) 0.0 ( 1 1 1) 0.75 ( 1 1 6) 10.00 TRANX 8.45 ( 1 1 4) 0.0 ( 1 1 1) 1.20383 ( 1 1 6) 7.042254 TRANY 8.45 ( 1 1 4) 0.0 ( 1 1 1) 1.20383 ( 1 1 6) 7.042254 TRANZ 1.49 ( 1 1 4) 0.0 ( 1 1 1) 0.1105796( 1 1 5) 13.45455 MULTX 1.0 ( 1 1 1) 1.0 ( 1 1 1) 1.0 ( 1 1 1) 1.0 MULTY 1.0 ( 1 1 1) 1.0 ( 1 1 1) 1.0 ( 1 1 1) 1.0 MULTZ 1.0 ( 1 1 1) 0.0 ( 1 1 9) 1.0 ( 1 1 1) 1.0 ROCKV 6028682. ( 1 1 6) 1939592. ( 1 4 4) 1939592. ( 1 4 4) 3.108222 MULTPV 1.0 ( 1 1 1) 1.0 ( 1 1 1) 1.0 ( 1 1 1) 1.0 MINPVV 0.1E-05 ( 1 1 1) 0.1E-05 ( 1 1 1) 0.1E-05 ( 2 1 1) 1.0 ------------------------------------------------------------------------------------------------------------------------------ ECLIPSE 300 Black-Oil Did you know that an ECLIPSE Black-oil license allows you to run ECLIPSE 300 in black-oil mode? ECLIPSE 300 sometimes gives you more information and insight into a problem than ECLIPSE 100. In-Place LGRs By default, ECLIPSE 100 will solve LGRs using local time stepping. With local time stepping, the global grid can take a large time step that is not limited by the size of the time step needed by the LGR. The algorithm used is shown below: The local time stepping method is semi-explicit and therefore is potentially unstable. The alternative is to use the LGRLOCK to “lock” the timesteps together and turn off the local time stepping, so that the global timesteps are the same as the local timesteps. With LGRLOCK, the method is fully implicit and unconditionally stable. The keyword LGRFREE will allow local timestepping to start again, if the convergence issue only occurs for a few timesteps. 6.3.2 SCAL data Another Story … Again, poor convergence was causing a model to run very slowly. Using the techniques described above, the problem was raced to one of the relative permeability curves. A plot of the original curve is shown below: The curve is smooth and looks reasonable. At first, the data in the relative permeability table also looked reasonable. The keyword was SOF3, and the data we were suspicious about was the krow as a function of Soil. -- SOIL SOF3 0.181 0.283 0.385 0.436 0.483 0.588 0.686 0.689 0.761 0.837 0.863 0.879 0.880 KROW 0 0.0001 0.0015 0.0124 0.0217 0.0939 0.3499 0.3501 0.7323 0.9887 0.9978 1 1 1 1 KROG 0 0.0001 0.0015 0.0124 0.0217 0.0939 0.3499 0.3501 0.7323 0.9887 0.9978 / The Soil and krow column see to be correct. The Soil column is monotonically increasing, and has values between 0 and 1. The krow column is also monotonically increasing, and has values between 0 and 1. We used the SCAL program to plot the derivative d(krow)/d(Soil). We expected that the derivatives would be smoothly varying, as krow was smoothly varying. We got an unexpected result: There are sharp changes in slope around Sw = 24% and around Sw = 31%. These value correspond to Soil = about 68% and 76%, so we checked the SOF3 table at these values. We found that at around 68% oil there are two oil saturations with very close values, 0.686 and 0.689, that have very similar krow values of 0.3499 and 0.3501. -- SOIL SOF3 0 0.283 0.385 0.436 0.483 0.588 0.686 0.689 0.761 0.837 0.863 0.879 0.880 KROW 0 0.0001 0.0015 0.0124 0.0217 0.0939 0.3499 0.3501 0.7323 0.9887 0.9978 1 1 1 1 KROG 0 0.0001 0.0015 0.0124 0.0217 0.0939 0.3499 0.3501 0.7323 0.9887 0.9978 / The SOF3 table does not need saturation values that close together, and in fact very closely spaced values in any table may cause ECLIPSE to do extra computing work as it has to find exactly which two point to interpolate between. In this case we can fix the problem by removing either one of the two entries: -- SOIL SOF3 0.181 0.283 0.385 0.436 0.483 0.588 0.689 0.761 0.837 0.863 KROW 0 0.0001 0.0015 0.0124 0.0217 0.0939 0.3501 0.7323 0.9887 0.9978 KROG 0 0.0001 0.0015 0.0124 0.0217 0.0939 0.3501 0.7323 0.9887 0.9978 0.879 0.880 1 1 1 1 / This single change was enough to produce smooth derivatives and to resolve the convergence problem. Removing very small kr values We can however also make an engineering improvement to the table. entry at Soil=0.283: -- SOIL SOF3 0.181 0.283 0.385 0.436 0.483 0.588 0.686 0.689 0.761 0.837 0.863 0.879 0.880 KROW 0 0.0001 0.0015 0.0124 0.0217 0.0939 0.3499 0.3501 0.7323 0.9887 0.9978 1 1 1 1 Consider the KROG 0 0.0001 0.0015 0.0124 0.0217 0.0939 0.3499 0.3501 0.7323 0.9887 0.9978 / That line in the table says that at 28.3% oil saturation the relative permeability of oil to water and to gas are both 1E-4. It is possible but unlikely that this 1E-4 is a value that has been measured experimentally. From a simulation point of view it makes more sense to set that value to 0. ECLIPSE has special code (known as the Appleyard Chop) that speeds up convergence at the point when a phase first becomes mobile. That code is wasted if the change in relative permeability is from 0 to 1E-4, but is very effective if the change is from 0 to a higher value. Critical saturations Consider a grid block with pore volume of 100 m3. Assume it contains 80 m3 of gas and 20 m3 of water. The water saturation is 20% Suppose 1. The pore volume decreases by 1% due to reservoir pressure changes 2. There is no flow in or out of that grid block. The gas will compress but water is relatively incompressible, so we will have approximately 79 m3 of gas and 20 m3 of water. The water saturation will increase from 20/100 = 20% to 20/99 = 20.2% As you produce/inject, the pressure in every block will change. Therefore the pore volume of each grid block will change. If the initial water saturation (SWL) in each grid block was 20%, Sw will now be >20% in many blocks. If SWCR=SWL=20% then all those blocks have mobile water even if Sw is only 20.0001%. This water WILL flow. Setting SWCR=SWL+0.01% will stabilize the model without changing the OOIP. The same applies to the critical saturation of any phase. The original SOF3 table is -- SOIL SOF3 0.181 0.283 0.385 0.436 0.483 0.588 0.689 0.761 0.837 0.863 0.879 0.880 KROW 0 0.0001 0.0015 0.0124 0.0217 0.0939 0.3501 0.7323 0.9887 0.9978 1 1 1 1 KROG 0 0.0001 0.0015 0.0124 0.0217 0.0939 0.3501 0.7323 0.9887 0.9978 / As discussed above, we could stabilize the oil by changing the 1E-4 values to zero. If for some reason we do not wish to do that, we could instead add a relative permeability of 0 at an oil saturation of 18.2% -- SOIL SOF3 0.181 0.182 0.283 0.385 0.436 0.483 0.588 0.689 0.761 0.837 0.863 0.879 0.880 KROW 0 0 0 0 0.0001 0.0015 0.0124 0.0217 0.0939 0.3501 0.7323 0.9887 0.9978 1 1 1 1 KROG 0.0001 0.0015 0.0124 0.0217 0.0939 0.3501 0.7323 0.9887 0.9978 / SCAL DATA: other issues Hysteresis • Does the model run better without hysteresis? • Are the drainage and imbibition end-points consistent? End-point scaling • Are the end-point values extreme (for instance SWL=>1)? • Are the end-points too close? • If you are using SCALECRS, is there a severe slope discontinuity at both ends of the 2-phase mobile regions? For example in an oil-water problem, do you have a discontinuity when (1-Sowcr) approaches SWU? In this case, try using SCALELIM to limit the value of (1-Sowcr) used in scaling. 6.3.3 PVT data Total Compressibility Check Common Warning (3) @--WARNING AT TIME 0.0 DAYS (19-OCT-1982): @ NEGATIVE COMPRESSIBILITY FOUND IN GAS @ PRESSURE TABLE 1 AND OIL PRESSURE @ TABLE 1 AT A SAMPLE PRESSURE VALUE @ 290.56207 . ADJUST SATURATED FLUID @ PROPERTY VALUES AT THIS PRESSURE. @ NEGATIVE COMPRESSIBILITIES OCCUR @ FOR GAS SATURATIONS LESS THAN 1.00000 Although this is only a warning, it is probably NOT safe to ignore the warning. In black oil models, both ECLIPSE 100 and ECLIPSE 300 check for positive compressibility of each single reservoir fluid as the PVT data is read (the formation volume factor must be a monotonically decreasing function of pressure assuming all other variables are held fixed). ECLIPSE also checks that mixtures of saturated oil and gas have a positive total compressibility even when there is mass transfer between the two phases. For example, increasing the pressure of a cell that contains both oil and gas will: • transfer some gas from the gas phase to the oil phase • swell the oil due to extra dissolved gas • decrease the remaining gas volume due to compression The oil volume will therefore increase with increasing pressure, and the gas volume will decrease with increasing pressure. The total (oil+gas) volume must however decrease. This will only happen if the reduction in the volume of gas is greater than the increase in the volume of oil. The total compressibility of the oil and gas system is given by the expression Ct = Sg*Ctg + So*Cto Ctg = [-dBg/dP + dRv/dP * (Bo-Rs*Bg) / (1-Rs*Rv)]/Bg Cto = [-dBo/dP + dRs/dP * (Bg-Rv*Bo) / (1-Rs*Rv)]/Bo Where Sg, So, Bg,Bo, Rs, and Rv have their usual meaning, and the pressure derivatives are taken along the saturation curve. In the special case of Rv=0, this simplifies to Ctg = [-dBg/dP ] / Bg Cto = [-dBo/dP + dRs/dP * Bg] / Bo For each PVT region in the simulation grid a pressure range is selected from the corresponding oil and gas PVT data tables to span the complete range of pressure data in the two tables. This pressure range is then subdivided into 30 equally spaced pressure nodes to evaluate the total hydrocarbon compressibility. At each pressure node two limiting compressibility values are calculated corresponding to Sg=0 and Sg=1: Ct = Ct,o at Sg = 0, So = 1 Ct = Ct,g at Sg = 1, So = 0 If either or both values are negative at a pressure node, you will get a warning message . The highest pressure used is the maximum of the oil pressure (maximum bubble point pressure in PVTO or maximum pressure specified in PVDO) and the gas pressure (maximum dew point pressure specified in PVTG or maximum gas pressure specified in PVDG). The maximum sample pressure at which total compressibility is checked can be overridden by choosing the second item in the PMAX keyword. If the range of sample pressures extends above the maximum bubble point specified in the PVTO table or above the maximum dew point specified in the PVTG table, then ECLIPSE is forced to extrapolate above the highest entered bubble point or dew point (see below). In this case, it is not unlikely that negative compressibilities could occur as a result of extrapolation. In ECLIPSE 100 only, if the 21st switch in the DEBUG keyword is set > 0, ECLIPSE 100 reports the total gas and oil compressibilities and as well as the saturated values of Rs, Bo and Rv,Bg for the sample range of pressures to the debug file (extension .DBG) in tabular form. Whenever any negative compressibilities are encountered in the data checking, this debug output is automatically switched on. You can use this information to diagnose any negative compressibilities reported by ECLIPSE 100. Typical output from DEBUG(21)>0 is: ============================================= Total compressibility report: Tcoil = total compressibility with Sg=0 Tcgas = total compressibility with So=0 Tc = Sg.Tcgas + So.Tcoil Minimum oil pressure 0.00000 Maximum oil pressure 8.23759 Minimum gas pressure 0.00000 Maximum gas pressure 196.1330 ============================================= Press 10.0000 16.7632 23.5264 30.2896 37.0528 43.8160 50.5792 57.3424 64.1057 70.8689 … Rssat 1.9000 2.1424 2.3148 2.4872 2.6596 2.8320 3.0044 3.1768 3.3493 3.5217 Rvsat 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Bg 1.0150 0.1336 0.0755 0.0493 0.0352 0.0272 0.0213 0.0155 0.0142 0.0138 Bo 1.0230 1.0395 1.0459 1.0523 1.0588 1.0654 1.0721 1.0400 1.0401 1.0401 DBg/DPo -1.0502 -0.0152 -0.0058 -0.0029 -0.0015 -0.0009 -0.0012 -0.0006 -0.0001 -0.0001 DBo/DPo 0.0254 0.0009 0.0009 0.0010 0.0010 0.0010 0.0010 0.0009 0.0009 0.0009 DRs/DPo 0.2039422 0.0254929 0.0254929 0.0254929 0.0254929 0.0254929 0.0254929 0.0254929 0.0254929 0.0254929 DRv/DPo 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 Tcoil 0.1774885 0.0023769 0.0009356 0.0002841 -0.0000670 -0.0002702 -0.0004215 -0.0005183 -0.0005506 -0.0005610 Tcgas 1.0346498 0.1135191 0.0772953 0.0589455 0.0421443 0.0343437 0.0543732 0.0397541 0.0045428 0.0044073 Since the total compressibility equation is: Ct = Sg*Ctg + So*Cto Then when P=70.8689 we can substitute the values of Ctg=Tcoil anf Ctg=Tcgas: Ct = Sg*0.0044073 – So*0.0005610 and Ct will be negative when So is close to 1 and Sg is close to 0. We have in this case of Rv=0, so the equations for Ctg and Cto are: Ctg = [-dBg/dP ] / Bg Cto = [-dBo/dP + dRs/dP * Bg] / Bo To correct the data in this case we can either decrease the absolute value of dBo/dP or increase the value of DRs/dP or of Bg. PVT DATA: Extrapolation Common Warning (4) @--WARNING AT TIME 0.0 DAYS (19-OCT-1982): @ NEW SOLUTION OBTAINED BY EXTRAPOLATION OF @ BOTH OIL AND GAS PVT TABLES This is only a warning, but it may NOT be safe to ignore it. In some cases, if the pressure is extrapolated severely, then the linearlyextrapolated formation volume factors may become unphysical. We recommend that the highest bubble point nodes in the PVTO table are constructed so avoid extrapolations above the highest Rs entered in the PVTO table during the simulation (similarly, in a run with vaporized oil present, it is recommended that the highest dew point nodes in the PVTG table are constructed to avoid extrapolations above the maximum Rv entered during the simulation). The keyword EXTRAPMS asks ECLIPSE 100 to warn you whenever extrapolations are made to PVT (or VFP) tables. In ECLIPSE 100, ALWAYS add keyword EXTRAPMS 4 / ECLIPSE 300 will warn you by default If insufficient PVT data is supplied, ECLIPSE may extrapolate the PVT table data to inaccurate or non-physical values. Question: If your initial reservoir pressure is 4800 psia, Is it OK to have PVT tables with a maximum pressure value of 4800 psia? Answer: No, because any injection is likely to be at a pressure greater than 4800 psia. If you are checking extrapolated values yourself, you should know that extrapolation is not linear in FVF or Viscosity. ECLIPSE stores PVT tables internally as the reciprocals of FVF and Viscosity* FVF. With EXTRAPMS value of 4, you may get extra information: @--WARNING @ @ AT TIME 0.0 DAYS (19-OCT-1982): NEW SOLUTION OBTAINED BY EXTRAPOLATION OF BOTH OIL AND GAS PVT TABLES Extrapolation in PVTO table number 1 Po-Pbub above highest value on undersaturated line Rs = 1.27027 MSCF/STB Po = 4820.658 PSIA Pbub = 4014.700 PSIA Cell (20,20,7) State 3 Sw = 0.12000 Sg = 0.00000 So = 0.88000 6.3.4 SOLUTION data Check that the initial solution is stable • If you are using EQUIL, it should be in equilibrium • For enumerated solutions, it depends on what you want to do • For RESTARTs, it won‟t be usually be in equilibrium To check for stability: 1. Run the model with no wells and no aquifers 2. Start with small timesteps and increase the timestep size 3. Check that there are no convergence problems 4. Check that field pressures and fluids-in-place do not change Plot SUMMARY quantities such as FPR, FOIP, … Fluids in place calculation There is a choice of methods to calculate the initial fluid in place values (the initial oil, water and gas saturations) in each grid block. The choice is selected using the ninth argument of the EQUIL keyword, which specifies the number of sub-layers (N) used to obtain the initial saturations: Center-point equilibration (N=0) The simulator sets the fluid saturations in each grid block according to the conditions at its center. This method gives a stable initial state since the phase pressure differences in the simulation are also taken between cell centers. It is however the least accurate method, particularly in cases where the fluid contact passes through large grid blocks. Horizontal fine grid equilibration (N<0) The top and bottom halves of each grid block are each divided into -N layers of equal thickness, and the saturations are determined locally in each layer. The phase saturations for the block are the average of the saturations in each layer. This option provides a more accurate calculation of the fluids in place, but may yield a solution that is not completely stable in the initial state. Tilted block fine grid equilibration (N>0) This option is similar to the previous option, but it takes into account the slope of each grid block. The top and bottom faces of the blocks are treated as planes that are tilted about their central points. The top and bottom halves of the grid block are each divided into N horizontal layers of equal thickness, but the thickness of the layers in the top half is generally different from the thickness of the layers in the bottom half. This is because the distance from the block center line to the highest corner is not the same as the distance from the block center line to the lowest corner, if the upper and lower faces have different tilts. The phase saturations in each block are calculated as a weighted average of the saturations in the layers, weighted according to the area of each layer that is enclosed within the block multiplied by the layer thickness. This option provides the most accurate calculation of fluids in place, but may yield a solution which is not completely stable in the initial state. With fine grid equilibration there is a redistribution of fluids between grid blocks near the contacts when the simulation begins, which occurs independently of any external driving force (wells etc.). The reason is that a steady state solution on the fine equilibration grid (in which each block is subdivided into several layers) is not necessarily a steady state solution on the coarser simulation grid. The redistribution of fluids can produce a significant transient when the simulation is started. In ECLIPSE 100 you can overcome this by setting the „quiescence switch' (the switch QUIESC in keyword EQLOPTS). If this switch is on, modifications are applied to the phase pressures to make the initial solution a true steady state. These pressure modifications are applied for the duration of the run. If EQUIL(9) is not zero, we may therefore have an unstable initial solution, so if the model with no wells and no aquifers is not stable: 1. Check the 9th argument of the EQUIL keyword 2. Make sure EQUIL(9) = 0 and run again. Alternatively in ECLIPSE 100 use QUIECENCE and run again. Still unstable? Are you using SWATINIT? • If the initial saturations defined in SWATINIT cannot be honoured then the model may not be stable • This can happen if you are using PPCWMAX with a large initial Sw a long way above the owc To check for aquifer stability: 1. Run the model with aquifers but with no wells 2. Use default values for the initial aquifer pressure 3. Start with small timesteps and increase the timestep size 4. Check that there are no convergence problems 5. Check that field pressures and fluids-in-place are stable Plot SUMMARY quantities such as FPR, FOIP, … 6. Do not worry about minor changes in FPR, FOIP, …. 6.3.5 SUMMARY section 1. Never use RPTONLY in the SUMMARY section if you are trying to find the cause of a problem. RPTONLY will only output summary information at report steps and not at timesteps. You may miss some important changes between timesteps. 2. Add PERFORMANCE to the SUMMARY section, and plot cpu time, elapsed time, timestep size, number of non-linears and number of linears. This can help you find the time during the simulation when the problems are happening, and may help you find the reason. 3. 6.3.6 SCHEDULE data Multi-Segment Wells Multi-Segment Wells are more difficult to solve than conventional wells. Add WSEGITER, with defaults If you have used TUNING (or any tuning keyword), make sure WSEGITER is after the TUNING keyword. WVFPEXP – If using VFP tables, can help for THP controlled wells especially with fluctuating gas production. WELSPECS(10) – turn off cross-flow. WELSPECS(12) – switching to AVG may reduce instability caused by explicit solution WELSEGS(8) in record 1 – if using drift flux model, try (a)switching to Homogeneous flow (b)(b) setting OPTIONS(77) WELSEGS(8) in record 2 – check if friction factor is too high. Lowering this may help convergence. OPTIONS(30) – governs methods aimed at improving convergence by chopping the increment in all the well‟s solution variables. OPTIONS(63) – can help if explicit calculation of the hydrostatic head between the well connections and the corresponding grid block centers causes oscillations Badly/incorrectly defined MSW definition can be the root cause of extreme convergence problems or even crashes! Defining MSWs by hand is extremely error-prone. Would only attempt for perfectly horizontal or vertical wells passing through regular box shaped grids. However, best to avoid doing by hand and use Petrel RE instead. There is currently no 3D visualisation for MSW wells, so don‟t be fooled by the Floviz display. There is no way to verify whether a manual definition is correct. 1. Don‟t do it by hand 2. If you do, don‟t rush it: Do it in stages: run the model with standard wells (no segments) to make sure simulation converges, then add the MSW definition. Draw the segments out on paper on the grid with exact numbers for depth and length before defining the MSWs Allow 1-2 days to define each well, more if it is complex. May be easier to do incremental definition than entering absolute numbers when defining the segments by hand. (WELSEGS item 5 of record 2) To have a segment per cell, need to have n segments, where n is the number of grid cells through which the well will pass. Segment 1 is the zero-tubing head segment. 7. Summary: Convergence ECLIPSE uses an iterative process based on Newton's method to solve the nonlinear equations Advance Timestep Linearize the Equations Iterate to solve the linear equations No Plug the linear solution into the non-linear equation Is the solution good? Yes “Non-linear iteration” The number of non-linear iterations is a guide to model convergence Non-linears per Timestep Guide 1 Very easy to converge 2 to 3 Easy to converge Increasingly difficult flow situation 4 to 9 > 10 Problem with model? Controls available from the Schedule section: TUNING sets timestepping, iteration and convergence criteria – TUNINGL is used for the LGRs in the model Guidelines: Timestepping controls need alteration fairly frequently Iteration controls seldom need adjustment Convergence controls need adjustment only in highly unusual circumstance 7.1 Convergence problem Checklist All sections Check all problem and warning messages. Does the model run better with no tuning keywords? Remove TUNING or TSCRIT or CVCRIT keywords and try to run. Is the model in equilibrium at the beginning of the simulation? Try running the model without any wells or any aquifers. Grid and EDIT sections Use MINPV, check LIMITS, check LGRs. Check whether you need to set a minimum pore volume Radial runs: Is the inner radius too small? Dual porosity runs: Is sigma too large LGRs: Are the cells too small? In ECLIPSE 100 try LGRLOCK PROPS section PVT data: Use EXTRAPMS in ECLIPSE 100. extrapolations In all simulators, avoid PVT SCAL data: Check the slopes of all relative permeability tables. SOLUTION section Is initial fluid distribution stable? SUMMARY section Remove RPTONLY Plot PERFORMANCE information SCHEDULE section Avoid VFP table extrapolations, use EXTRAPMS in ECLIPSE 100. Take care when using MSWs
© Copyright 2024