水资源短缺风险综合评价.docx

上传人:小飞机 文档编号:1809040 上传时间:2022-12-20 格式:DOCX 页数:31 大小:321.55KB
返回 下载 相关 举报
水资源短缺风险综合评价.docx_第1页
第1页 / 共31页
水资源短缺风险综合评价.docx_第2页
第2页 / 共31页
水资源短缺风险综合评价.docx_第3页
第3页 / 共31页
水资源短缺风险综合评价.docx_第4页
第4页 / 共31页
水资源短缺风险综合评价.docx_第5页
第5页 / 共31页
点击查看更多>>
资源描述

《水资源短缺风险综合评价.docx》由会员分享,可在线阅读,更多相关《水资源短缺风险综合评价.docx(31页珍藏版)》请在三一办公上搜索。

1、2013高教社杯全国大学生数学建模竞赛承 诺 书我们仔细阅读了中国大学生数学建模竞赛的竞赛规则.我们完全明白,在竞赛开始后参赛队员不能以任何方式(包括电话、电子邮件、网上咨询等)与队外的任何人(包括指导教师)研究、讨论与赛题有关的问题。我们知道,抄袭别人的成果是违反竞赛规则的, 如果引用别人的成果或其他公开的资料(包括网上查到的资料),必须按照规定的参考文献的表述方式在正文引用处和参考文献中明确列出。我们郑重承诺,严格遵守竞赛规则,以保证竞赛的公正、公平性。如有违反竞赛规则的行为,我们将受到严肃处理。我们授权全国大学生数学建模竞赛组委会,可将我们的论文以任何形式进行公开展示(包括进行网上公示,

2、在书籍、期刊和其他媒体进行正式或非正式发表等)。我们参赛选择的题号是(从A/B/C/D中选择一项填写): 我们的参赛报名号为(如果赛区设置报名号的话): 所属学校(请填写完整的全名): 参赛队员 (打印并签名) :1. 2. 3. 指导教师或指导教师组负责人 (打印并签名): 日期: 年 月 日赛区评阅编号(由赛区组委会评阅前进行编号):2013高教社杯全国大学生数学建模竞赛编 号 专 用 页赛区评阅编号(由赛区组委会评阅前进行编号):赛区评阅记录(可供赛区评阅时使用):评阅人评分备注全国统一编号(由赛区组委会送交全国前编号):全国评阅编号(由全国组委会评阅前进行编号):水资源短缺风险综合评价

3、摘要水资源短缺问题是影响我国发展的重大问题,本文针对水资源短缺风险问题找出了主要风险因子,建立了水资源短缺风险评价模型,对水资源短缺风险进行等级划分,并提出相应的有效措施规避风险。对于问题一,我们建立主成分和灰色关联度分析模型,分析附表和相关资料,先确立了北京市水资源短缺风险的风险因素主要包括自然因素,即降雨量和常住人口,和社会因素,即农业用水,工业用水,第三产业及生活其他用水,污水处理率,城市绿化覆盖率。然后利用主成分分析得到个各个因子的贡献率,再利用灰色关联度分析,得到各个因子与缺水量的关联度的大小,基本与主成分分析一致,最后得到主要风险因子。对于问题二,我们用综合评价的模型,对模型一所确

4、定的主要风险因子做相应高优和低优指标处理,并对北京市水资源短缺进行风险等级划分。最后对主要风险因子进行调控,来降低风险等级。对于问题三,我们建立 模型,要对北京市未来两年水资源的短缺风险进行预测,我们通过对主要风险因子进行预测,并对预测模型进行后验差检验,然后再用模型,给未来的两年划分风险等级。对于问题四,我们通过分析上面的数据和查找相关资料,给北京市水行政主管部门写一份建议报告。关键词:主成分分析 灰色关联度分析 模型 模型 后验差检验 一、问题重述新中国成立以来,党和政府领导全国人民进行了大规模的水利建设,取得了巨大成就,但水的问题仍然困扰着我国经济社会的发展。1998年中国水资源公报中的

5、数据表明:中国面临的水资源形势依然严峻。水的问题已经成为制约我国绝大部分地区经济社会可持续发展的重要因素。水资源,是指可供人类直接利用,能够不断更新的天然水体。主要包括陆地上的地表水和地下水。风险,是指某一特定危险情况发生的可能性和后果的组合。水资源短缺风险,泛指在特定的时空环境条件下,由于来水和用水两方面存在不确定性,使区域水资源系统发生供水短缺的可能性以及由此可能产生相关的的损失。北京是世界上水资源严重缺乏的大都市之一,目前,北京水资源年均用水总量达36亿立方米,而年均水资源总量仅有21亿立方米,缺口达15亿立方米。北京的年均水资源量平摊到每个人身上,甚至不足100立方米。为全国人均的1/

6、8,世界人均的1/30,属重度缺水地区,为了维持到现在,北京用尽了各种方法,除了寻找多种水源“开源”,还通过调整用水结构、再生水利用等办法“节流”。用三句话概括,就是农业用水负增长,工业用水零增长,生活用水适度增长。附表中所列的数据给出了1979年至2000年北京市水资源短缺的状况。北京市水资源短缺已经成为影响和制约首都社会和经济发展的主要因素。政府采取了一系列措施, 如南水北调工程建设, 建立污水处理厂,产业结构调整等。但是,气候变化和经济社会不断发展,水资源短缺风险始终存在。如何对水资源风险的主要因子进行识别,对风险造成的危害等级进行划分,对不同风险因子采取相应的有效措施规避风险或减少其造

7、成的危害,这对社会经济的稳定、可持续发展战略的实施具有重要的意义。北京2009统计年鉴及市政统计资料提供了北京市水资源的有关信息。利用这些资料和我们可获得的其他资料,讨论了以下问题:问题一:影响水资源的因素很多,例如:气候条件、水利工程设施、工业污染、农业用水、管理制度,人口规模等。我们评价判定北京市水资源短缺风险的主要风险因子是什么。问题二:建立一个数学模型对北京市水资源短缺风险进行综合评价, 作出风险等级划分并陈述理由。并通过相关数据对主要风险因子,如何进行调控,使得风险降低。问题三:建立预测模型,对北京市未来两年水资源的短缺风险进行预测,并提出应对措施。问题四:以北京市水行政主管部门为报

8、告对象,写一份建议报告。二 、基本假设(1)用水量是农业用水,工业用水,第三产业及生活其他用水的总和。模型的分析与建立均以北京市为基础。(2)由于水资源的不稳定性,在计算中排除某一差异较大的数,对模型无影响。模型的分析与建立均以北京市为基础。(3)引起水资源短缺的因素除我们选择的7个因子外,其他的因子影响都很小,可以忽略不计。(4)在一年时间内,各风险因子对北京市水资源短缺的影响是一定的。(5)模型的分析与建立均以北京市为基础。三、符号说明:分辨率;:矩阵相关系数;: 评价对象;: 评价指标排成;: 各组数据的秩;: 相应权重;: 称为发展灰数;: 称为内生控制灰数四、问题分析对于问题一,我们

9、建立的是用主成分和灰色关联度分析法对水资源短缺的风险因子进行综合评价。通过收集资料和分析数据,得到影响水资源短缺风险因素,主要有社会因素和自然因素,社会因素又可以细分为农业用水,工业用水,第三产业及生活其他用水,污水处理率,城市绿化覆盖率。自然因素可以细分为降雨量和常住人口。主成分分析中得到每个因子在各个主成分中的贡献率,根据贡献率的大小提取主要因子。然后再结合灰色关联度分析来求得各个因子与缺水量的关联度,对主成分分析中得到的主要风险因子进行检验。对于问题二,基于问题一所确立的主要风险因子,和对低优指标相应处理,利用模型对北京市水资源短缺风险进行综合评价。根据对应的概率单位值进行相应等级划分。

10、通过控制影响北京市水资源短缺的主要风险因子的权重,来改变风险等级。对于问题三,我们用 模型,分别预测了未来两年主要风险因子,即污水处理率,农业用水,第三产业及生活其他用水,降雨量,常住人口,在利用问题二的模型,得到未来两年的风险等级。对于问题四,通过对水资源短缺风险因子的分析,和相关的预测,我们可以采取一些措施,降低水资源短缺的风险等级。五、模型的建立与求解5.1问题一模型的建立与求解通过查找相关资料,得到影响水资源短缺风险因素,主要有社会因素即为农业用水,工业用水,第三产业及生活其他用水,污水处理率,城市绿化覆盖率。和自然因素即降雨量和常住人口。如图1所示。水资源短缺风险因素社会因素自然因素

11、农业用水污水处理率常住人口城市绿化覆盖率降雨量第三 产业及其他用水图1风险因子分布图工业用水5.1.1主成分分析的研究原理在数据分析工作中,常常需要把很复杂的数据集简化,即将 个指标所构成的 维简化为一位系统。主成分分析就是多个指标化为少数几个综合指标的一种统计分析方法。在多指标(变量)研究中,往往由于变量个数太多,且彼此之间存在着一定的相关性,因而使得所观测的数据在一定程度上有信息的重叠。当变量较多时,在高维空间中研究样本的分布规律就更复杂。主成分分析采取一种降维的方法,找出几个综合因子尽可能地反映原来变量的信息量,而且彼此之间互不相关,从而达到简化的目的。主成分分析的计算步骤如下:设观测样

12、本矩阵为(为样本数, 为变量数):为使该样本集在降维中所引起的平方误差最小,必须进行两方面的工作:一是进行坐标变换,即用雅可比方法求解正交变换矩阵;二是选取个主成分.将原始数据进行标准化处理,即对样本集中元素:作变换 (i=1,2,n;k=1,2, )作变换,即,其中:为参与评价的指标个数。主成分分析的明显特征是每个主分量依赖于测量初始变量所用的尺度,当尺度改变时,会得到不同的特征值。克服这个困难的方法是对初始变量进行以上标准化处理,使其方差为1.计算样本矩阵的相关系数矩阵。对应相关系数矩阵,计算特征值以及各个主成分的贡献率。再把贡献率大的因子找出来,从而起到筛选的作用。5.1.2 主成分分析

13、法求解过程在附录1中给出了1979到2010年北京水资源短缺风险评价的原始数据,我们利用 将这些数据进行标准差标准化处理后计算其相关系数矩阵(见表1) 表1 相关系数矩阵污水处理率(%)农业用水(亿立方米)工业用水(亿立方米)第三产业及生活其他用水(亿立方米)降雨量(毫米)城市绿化覆盖率(%)常住人口(万人)污水处理率(%)1.000-0.648-0.9240.877-0.3800.8670.946农业用水(亿立方米)-0.6481.0000.506-0.710-0.019-0.712-0.722工业用水(亿立方米-0.9240.5061.000-0.8130.397-0.814-0.879第

14、三产业及生活其他用水(亿立方米0.877-0.710-0.8131.000-0.2920.9640.958降雨量(毫米)-0.380-0.0190.397-0.2921.000-0.318-0.300城市绿化覆盖率(%)0.867-0.712-0.8140.964-0.3181.0000.946常住人口(万人)0.946-0.722-0.8790.958-0.3000.9461.000再由相关系数矩阵计算特征值和各个主成分的贡献率和累积贡献率如表2,表3所示。 表2 特征值和贡献率和累积贡献率成分初始特征值提取平方和载入合计方差的%累积 %合计方差的%累积 %15.24374.90074.90

15、05.24374.90074.90021.05014.99489.8941.05014.99489.89430.3815.44095.3340.3815.44095.33440.2072.95898.2920.2072.95898.29250.0650.93399.2240.0650.93399.22460.0340.48299.70670.0210.294100.000 表3 成分矩阵12345污水处理率(%).958-.067.134-.167.157农业用水(亿立方米)-.742-.486.417.194.042工业用水(亿立方米)-.905.164-.308.186.155第三产业及生

16、活其他用水(亿立方米).962.075-.008.228.021降雨量(毫米)-.374.878.298.023.004城市绿化覆盖率(%).960.054-.039.228-.085常住人口(万人).985.058.060.044.085由表2、3可以看出,第一、第二主成分的累计贡献率已达到了89.894%,我们在第一和第二主成分中选择因子,各个因子系数的大小可以反映因子对主成分的贡献率,所以我们选择了农业用水,第三产业及生活其他用水,降雨量,城市绿化覆盖率,常住人口,作为水资源短缺的主要风险因子。 5.1.3 灰色关联度的原理和求解灰色关联度分析是分析系统中各因素关联的程度的方法,在计算关

17、联度之前,需先计算关联系数。关联系数:则关联系数定义为: 其中: 为第个点和 的绝对误差;为两级最小差;为两级最大差; 称为分辨率,0 1,若 越小,关联系数间差异越大,区分能力越强。一般取 =0.5。对单位不一,初值不同的序列,在计算相关系数前应首先进行初始化,即将该序列所有的数据分别除以第一个数据。利用上述方法,我们缺水量来衡量水资源短缺程度,所以把缺水量作为参考序列,各个因子作为比较序列,得到各个因子与缺水量的关联度(如表4所示)。表4 各因子关联度值因子污水处理率农业用水工业供水第三产业及生活其他用水降雨量城市绿化率常住人口0.55010.65200.58920.69370.65020

18、.69910.7042根据关联度的大小,我们同样可以选择农业用水,第三产业及生活其他用水,降雨量,城市绿化率,常住人口。5.1.4 结果说明用主成分分析评价影响水资源短缺的主要风险因子和用灰色关联度分析方法对影响水资源短缺的主要风险因子,其结果基本一致,说明对于北京水资源短缺的主要风险因子为农业用水,第三产业及生活其他用水,降雨量,城市绿化率,常住人口。5.2问题一模型的建立与求解5.2.1 原理秩和比()指的是表中行(或列)秩次的平均值,是一个非参数计量的综合指数,具有01区间连续变量的特征,其所有比较组秩和比之和为 ;如果编秩不按照经典的秩变换方法,各组秩和比的合计可能不为 。其基本思想是

19、在一个行(评价对象)列(个评价指标或等级)矩阵中,通过秩转换,获得无量纲的统计量,以值对评价对象的优劣进行排序,进而根据比较组数的多少,进行分档处理(比组数较多)或进行 平方根反正弦变换值可信区间处理(比较组数较少)。5.2.2 模型的求解(1)确定评价对象和综合评价因子和编秩将 个评价对象和个评价指标排成行列的原始数据矩阵。高优指标从小到大编秩,低优指标从大到小编秩,相同者取平均值。本题我们确定5个主要影响因素:农业用水、第三产业及其他用水、常住人口、污水处理率、降雨量。从附录表1中提取它们对应的32年原始数据。(2)计算和通过软件计算各组数据的秩()和相应权重(),公式如下: (1) (2

20、)应用上一步的5组原始数据,通过公式(1)、(2)得到如下结果:= 0.3452 0.5080 0.4369 0.4156 0.4655 0.3901 0.2009 0.4369 0.3030 0.4831 0.5742 0.4972 0.5458 0.5813 0.5833 0.5334 0.5618 0.4998 0.5314 0.4829 0.6199 0.5959 0.6061 0.5583 0.5910 0.5731 0.5997 0.6211 0.5906 0.5675 0.6087 0.5914 (3)的分布的分布是指用概率单位 表达的定值的向下累计频率 。编制频数分布表 ,列出

21、各组频数、累计频数 ;确定各组的秩次 及平均秩次 计算向下累计频率 ,修正最后一个累计频率,最后一个累积频率按 计算;将百分率 换算为概率单位 ,它 为百分率 对应的标准正态离差u加5。表5 的分布累计频数0.2091110.03133.13730.3001220.06253.46590.3421330.09383.68200.3911440.12503.84970.4161550.15633.99000.4391660.18754.11290.4391770.21884.22360.4651880.25004.32550.4891990.28134.42090.481110100.31254

22、.51120.492111110.34384.59770.498112120.37504.68140.500113130.40634.76280.534114140.43754.84270.534115150.46884.92160.548116160.50005.00000.553117170.53135.07840.568118180.56255.15730.565119190.59385.23720.571120200.62505.31860.572121210.65635.40230.583122220.68755.48880.583123230.71885.57910.5961242

23、40.75005.67450.590125250.78155.77640.594126260.81255.88710.599127270.84386.01000.597128280.87506.15030.601129290.90636.31800.607130300.93756.53410.619131310.96886.86270.621132320.99227.4176(4)计算回归方程以累计频率对应的概率单位值为自变量 ,值为变量 ,计算回归方程。根据表5,利用软件计算回归方程结果为方程拟合优度为0.9147,说明通过检验。(5)分档排序根据表5常用分档数与百分位数及概率单位Y,以值对

24、评价对象进行分档排序表6 年常用分档数与百分位数及概率单位挡数百分位数挡数百分位数315.866以下710.02715.86633.3684.13467.00346.681以下89.9736.68198.3525081.222以下93.3191.22253.593以下6.6813.59322.66327.4255072.57577.33796.40793.31962.275以下98.6782.27590.990以下15.8660.99504.74684.13415.86697.72537.0771.168以下62.931.16884.134根据表5、6,结合本题需要,我们确定采用5档的分法 ,

25、结果如表7所示:表7 结果分析等级等级对应年份对应年份总数的范围一级198516.80(6)秩和检验秩和检验用于检验分档是否合理,同等档级之间会有不和谐现象,用秩和检验可以排除此类问题的干扰。 秩和检验值,当 ,说明分档和谐,否则不和谐。我们利用程序(见程序4)计算出各等级秩和检验值如下表表8 各等级对应秩和检验值等级检验值是否和谐一级1和谐二级0.9825和谐三级0.9199和谐四级0.7809和谐五级0.9779和谐通过表8我们可以知道各等级是和谐的,因此我们的风险等级划分合理。5.3. 问题三的模型建立与求解5.3.1模型原理在对模型对实际问题进行预测时发现,影响预测精度的因素有很多时,

26、其中发展灰数增长率较小时,预测精度越高,当原始数据增长率变化较大时,即模型预测方法中地发展灰数较大时,预测精度低。模型预测,具体步骤如下:根据原始数据,建立时间序列:通过将原始数据进行一次性累加生成新序列:其中:建立模型相应的微分方程为:构造累加数据矩阵和常数向量用最小二乘法求得灰色参数其中:称为发展灰数,称为内生控制灰数。设为待估参数向量,解的:求解微分方程,即可得预测模型:5.3.2模型的求解通过建立的 模型,因为降雨量存在着随机性,所以我们对农业用水,第三产业及生活其他用水,城市绿化率,常住人口,进行预测拟合,利用编程(见附录程序3),画出了拟合散点图,并得到了拟合方程。 圆圈点为原始值

27、,菱形点为预测值。图2 农业用水的拟合图从图中可以看到原始值和预测值的都是以递减的方向变动,原始值在预测值的上下均匀分布。对农业用水的拟合方程: 图3 第三产业及生活其他用水的拟合图从图中可以看到原始值和预测值的都是以递增的方向变动,原始值在预测值的上下均匀分布。对第三产业及生活其他用水的 GM (1,1) 拟合方程图4 城市绿化率的拟合图从图中可以看到原始值和预测值的都是以递增的方向变动,原始值在预测值的上下均匀分布。对城市绿化率的 GM (1,1)拟合方程:图5 常住人口的拟合图从图中可以看到原始值和预测值的都是以递增的方向变动,原始值在预测值的上下均匀分布。对常住人口的拟合方程:5.3.

28、3模型的检验对模型进行后验差检验(1)计算原始序列标准差:(2)计算绝对误差序列的标准差:(3)计算后验差比值:(4)计算小误差概率: 等级0.950.800.700.65勉强合格=0.65不合格表9 检验值对照表用上述检验方法对主要风险因子进行检验,结果如表10所示表10 风险因子检验风险因子农业用水0.65280.8438勉强合格第三产业及生活其他用水0.23431好城市绿化率0.17841好常住人口0.18810.9688好主要风险因子基本上通过了后验差检验,所以可以用拟合方程来预测未来两年的主要风险的值,结果如表11所示。表11 2011年和2012年风险因子预测值风险因子农业用水第三

29、产业及其它用水常住人口污水处理率201111.321.09411848103.5201210.9922.08101893115.1 再根据的MATLAB程序(如附录程序2)计算出2011年和2012年的,查附件表2百分数与概率单位对应值可以得到如下结论:表12 2011和2012风险等级预测年份风险等级20110.57075.1789第三等级20120.60695.2715第三等级六、模型的评价质合比法因为数据太少,主观因素太大,可能会造成误差。而且质合检验缺点是相同秩次较多时,统计量要经过多次校正。关联度分析的缺陷在于基本思想欠缺准确和规范内涵欠全面两个方面,而本文也无法避开这两个方面,这一

30、点有待改进。 方法本身是一种确定的外推,在处理数据、拟合曲线和得到的拟合曲线的过程中都不考虑随机误差。参考文献1 苏金明,王永利,MATLAB7.0使用指南M,北京:电子工业出版社,2004.2 江世宏,MATLAB语言与数学实验M,北京:科学出版社,2007.3 赵静,但琦,数学建模与数学实验M,高等教育出版社,2008.1.4 冉启康,张振宇,张立柱,常用数学软件教程,人民邮电出版社,2008.10.5 若毛吉,应用主成分和灰色关联度分析法评价共和盆地天然草地产草量的影响,6 吴清平,张丹,秩和比法和几种常用评价方法在医疗质量评价中的比较J, 中国医院统计,2003.03.25.7 冯长春

31、, 侯玉亭,城镇土地评价中主成分分析法的应用J, 中国国土资源经济, 2007.07.15.8 何国华,区域物流需求预测及灰色预测模型的应用, 北京交通大学学报J, 2008.01.05.9 田志红,我国水资源保护存在的问题与对策研究J,科技情报开发与经济,2010年02期.附录表1 1979到2010的相关数据年份污水处理率(%)农业用水(亿立方米)工业用水(亿立方米)第三产业及生活他用水(亿立方米)降雨量(毫米)197910.224.1814.374.37718.419809.431.8313.774.94380.7198110.831.612.214.3393.2198210.928.8

32、113.894.52544.4198310.231.611.244.72489.919841021.8414.3764.017488.819851010.1217.24.3972119868.919.469.917.18665.319877.79.6814.017.26683.919887.421.9914.046.4673.319896.624.4213.776.45442.219907.321.7412.347.04697.319916.622.711.97.43747.919921.219.9415.5110.98541.519933.120.3515.289.59506.719949.

33、620.9314.5710.37813.2199519.419.3313.7811.77572.5199621.218.9511.769.3700.919972218.1211.111.1430.9199822.517.3910.8412.2731.719992518.4510.5612.7266.9200039.416.4910.5213.39371.120014217.49.212.3338.920024515.57.511.6370.4200350.113.88.413.6444.9200453.913.57.713.4483.5200562.413.26.814.5410.720067

34、3.812.86.215.3318200776.212.45.816.6483.9200878.9125.217.9626.3200980.3125.218.3480.620108111.45.118.8522.5表2 百分数与概率单位对照表%0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 0 -1.9098 2.1218 2.2522 2.3479 2.4242 2.4879 2.5427 2.5911 2.6344 1 2.6737 2.7096 2.7429 2.7738 2.8027 2.8299 2.8556 2.8799 2.9034 2.9251

35、2 2.9463 2.9665 2.9859 3.0046 3.0226 3.0400 3.0569 3.0732 3.0890 3.1043 3 3.1192 3.1337 3.1478 3.1616 3.1759 3.1881 3.2009 3.2134 3.2256 3.2376 4 3.2493 3.2608 3.2721 3.2831 3.2940 3.3046 3.3151 3.3253 3.3354 3.3454 5 3.3551 3.3648 3.3742 3.3836 3.3928 3.4018 3.4107 3.4195 3.4282 3.4268 6 3.4452 3.4

36、536 3.4618 3.4699 3.4780 3.4859 3.4937 3.5015 3.5091 3.5167 7 3.5242 3.5316 3.5389 3.5462 3.5534 3.5606 3.5675 3.5745 3.5813 3.5882 8 3.5949 3.6016 3.6083 3.6148 3.6213 3.6278 3.6342 3.6405 3.6468 3.6531 9 3.6592 3.6654 3.6715 3.6775 3.6835 3.6894 3.6953 3.7012 3.7070 3.7127 10 3.7184 3.7241 3.7298

37、3.7354 3.7409 3.7464 3.7519 3.7547 3.7625 3.7681 11 3.7735 3.7788 3.7840 3.7893 3.7945 3.7996 3.8048 3.8099 3.8150 3.8200 12 3.8250 3.8300 3.8350 3.8399 3.8448 3.8497 3.8545 3.8593 3.8641 3.8689 13 3.8736 3.8783 3.8830 3.8877 3.8923 3.8969 3.9015 3.9061 3.9107 3.9152 14 3.9197 3.9242 3.9268 3.9331 3.9375 3.9419 3.9463 3.9

展开阅读全文
相关资源
猜你喜欢
相关搜索

当前位置:首页 > 生活休闲 > 在线阅读


备案号:宁ICP备20000045号-2

经营许可证:宁B2-20210002

宁公网安备 64010402000987号